Drug-target affinity prediction method based on multi-scale mixed attention network

Abstract

The invention relates to the field of drug-target affinity prediction, in particular to a drug-target affinity prediction method based on a multi-scale mixed attention network. According to the method, global residue characteristics can be learned by protein pocket residues by utilizing a self-attention mechanism, so that the global characteristics of the protein are learned. In addition, a cross-modal feature fusion mechanism is utilized, protein and drug structural features are enhanced in sequence features, cross-modal feature fusion in entities is enhanced, meanwhile, a cross-entity interaction module is utilized to recognize key atoms of drug molecules or key residues of protein, and therefore the overall performance of the model is improved. Finally, multi-level protein feature extraction is further introduced to further extract protein features.

Description

Technical Field

[0001] The present invention relates to the field of drug-target affinity prediction, and in particular to a drug-target affinity prediction method based on a multi-scale hybrid attention network. Background Art

[0002] Discovering drugs that target specific proteins plays a key role in modern drug research. Among the many evaluation indicators, drug-target affinity (DTA) is a commonly used indicator, which reflects the binding strength between drug molecules and protein targets, thus facilitating drug screening targeting proteins. However, the process of determining DTA by traditional experimental measurement methods is usually time-consuming and costly. Despite great efforts in energy-based computational methods such as molecular dynamics simulations and molecular docking, these tools still have problems in balancing efficiency and accuracy. Therefore, the development of an efficient and effective computer analysis tool remains the top priority for current drug discovery.

[0003] In recent years, deep learning-based methods have been widely used in drug-target affinity (DTA) tasks. These methods are different from traditional molecular dynamics simulations in that they predict affinity by learning from known samples, bypassing complex energy rules. According to the availability of binding conformations of drug-target complexes, existing methods are divided into two categories: non-interaction methods and interaction-based methods. Non-interaction methods take drugs and proteins as inputs separately, but fail to combine specific conformational information for analysis; interaction-based methods take specific local binding regions of the complex as input. For non-interaction methods, drugs and proteins are represented separately due to the lack of specific binding conformations. In general, these inputs usually use one-dimensional descriptors, such as SMILES (Simplified molecular input line entry specification) for drugs and amino acid residue sequences for proteins; or two-dimensional descriptors, such as molecular graphs for drugs and residue contact graphs for proteins. Then, features are extracted using convolutional neural networks (CNNs) or graph neural networks (GNNs) and their countless variants. At present, various methods use deep learning technology to predict the interaction between drugs and target proteins. These methods can be roughly divided into two categories: one is the non-interaction method, which takes drugs and proteins as inputs respectively, and focuses on predicting affinity by learning known samples, usually using one-dimensional or two-dimensional descriptors, such as SMILES and molecular graphs of drugs. However, they cannot capture fine-grained structural information, such as specific binding sites and physical forces between molecules. The other is the interaction-based method, which pays more attention to extracting structural features from the binding conformation of the drug-target complex, especially the information of the local binding area, and improves feature extraction through technologies such as three-dimensional CNN and GNN. These methods can provide more detailed structural descriptions, but the computational complexity is high. In general, the non-interaction method is lightweight and efficient, suitable for large-scale screening; while the interaction-based method provides more accurate structural information and is suitable for fine analysis. Although the current interaction-based methods have achieved good performance, there are still some potential disadvantages: First, there is a lack of global awareness. Most current methods focus on mining the features of local binding regions (including pocket regions and drug molecules), while ignoring the potential association between pocket residues and distal residues. Second, insufficient representation. For most interaction-based methods, structural information is always fully utilized due to the availability of drug-target complexes; however, features of sequential patterns, including residue sequences of proteins and SMILES of drugs, are either ignored or simply concatenated. Third, the interaction is insufficient, and fine-grained interactions between drugs and targets are missing. The representations of drugs and targets are learned in an isolated manner, and the learned features are simply concatenated for the final prediction.

[0004] The present invention proposes a multiscale hybrid attention network for drug-target affinity prediction, referred to as MHAN-DTA (Multiscale HybridAttention Network for the DTA prediction task), which aims to solve the problem of insufficient feature mining and thus improve the prediction performance. In order to give the model global perception ability, a pocket-oriented feature aggregation and extraction module based on the self-attention mechanism is developed, and a hierarchical strategy is applied to the target protein. The present invention further introduces a cross-modal fusion module and a cross-entity interaction module to mine multi-scale intramolecular and intermolecular features within the binding site. Summary of the invention

[0005] Drug-target affinity prediction is a key challenge in the drug discovery process. Recent advances have demonstrated the great potential of deep learning for affinity prediction; however, existing methods fail to learn sufficient representations of drug-target complexes, resulting in poor performance.

[0006] Drug-target affinity prediction refers to inputting the features of proteins and drugs into the deep learning model, so that the model can quickly predict the affinity between drugs and protein targets, significantly reducing the economic and time costs of traditional biochemical experiments. In view of the shortcomings of previous drug-target affinity prediction models, the technical problems solved by the present invention are as follows:

[0007] (1) A pocket-oriented feature aggregation and extraction module was developed to model the association between local pocket residues and remote global residues, thus improving the global perception of protein representation. At the same time, a hierarchical strategy was adopted to capture the global representation of the target protein.

[0008] (2) A cross-modal feature fusion module is developed, in which the sequence modalities of proteins and drugs are fused through a cross-attention mechanism, making full use of intramolecular features to generate sufficient representations.

[0009] (3) A cross-entity feature interaction module is developed, where each dimension of the learned drug or pocket representation interacts with each feature dimension of another entity. This thorough interaction further improves the intermolecular features.

[0010] By collecting the sequence and structural information of proteins and drugs, this paper proposes a cross-modal fusion and cross-entity interaction mechanism based on a multi-scale hybrid attention network. This method effectively integrates the sequence information, structural features and interaction relationship between proteins and drugs during the training process of the drug-target affinity prediction model.

[0011] In order to make up for the shortcomings of the prior art, the present invention provides a drug-target affinity prediction method based on a multi-scale hybrid attention network.

[0012] The present invention is implemented by the following technical scheme: a drug-target affinity prediction method based on a multi-scale hybrid attention network, comprising the following steps:

[0013] Step 1: Construction of drug SMILES and molecular graphs as well as residue sequence and pocket graphs of target proteins

[0014] Drug molecules are usually represented by their natural sequence descriptors SMILES, which can be defined as Seq d ={c 1 ,c 2 ,...,c m}, where each element c i is the character that represents the chemical bond or atom in the SMILES string. i There are 63 characters in total, so Seq d It is encoded as a numerical sequence consisting of numbers in [0,63). After obtaining SMILES, the molecular graph is calculated through SMILES. First, the drug SMILES is converted into a molecular graph using RDKit. RDKit is an open source toolkit for chemical informatics. It is based on the operation of 2D and 3D molecules of compounds and uses machine learning methods to generate compound descriptors, fingerprints, calculate compound structural similarity, 2D and 3D molecular displays, etc. It can be called and used based on the Python language; the drug graph descriptor Gra d = {V d ,E d},V d represents the set of nodes containing atoms in a molecule, E d is the key edge set. Therefore, the initial feature matrix can be obtained and the adjacency matrix Where N d =|V d | represents the number of atomic nodes in the graph, D d =63 represents the characteristic dimension of the atom, represents the real number domain, that is, all elements of this matrix are real numbers. All indicate this meaning.

[0015] Proteins are naturally composed of sequences of amino acid residues. Therefore, the target protein can be directly represented as a residue sequence Seq p = {r 1 ,r 2 ,...,r n}, where each element r i Represents a residue. There are 20 types of residues in total, plus one filler value type, so Seq pis encoded into a digital sequence consisting of numbers in [0,21). For graph construction, in order to better focus on the local binding mode, the present invention constructs a graph for a local protein pocket rather than the entire protein. Residues whose distance to the nearest drug molecule atom is less than a certain threshold δ are considered pocket residues. The threshold δ hyperparameter is set to 5 angstroms in the present invention and can be further adjusted to 3 to 10 angstroms based on experimental measurements. Then, an edge set E is generated based on the Euclidean distance between the Cα atom pairs of any two residues. pocket Only when the distance is less than a certain threshold value θ, two residual nodes will be connected. The threshold value θ is a hyperparameter, which is set to 6 angstroms in the present invention and can be further adjusted to 3 to 10 angstroms according to experimental measurements. Based on the above process, the descriptor Gra of the protein pocket can be generated. pocket = {V pocket ,E pocket},V d represents the node set containing the residues in the protein pocket, E d Represents the edge set of the protein pocket residue contact graph. And the adjacency matrix Where N pocket =|V pocket | indicates the number of pocket residues.

[0016] Step 2: Use multi-head self-attention mechanism to achieve global feature aggregation

[0017] For the protein branch, the residue sequence Seq obtained in step 1 is p Input embedding layer to generate initial representation Where n is the length of the protein sequence and e = 128 is the embedding dimension. Next, the generated data is encoded with the same dimension through the sinusoidal encoding scheme, i.e. Here is the encoding process:

[0018] PE (i,2j) = sin(i / 10000 2j / e )

[0019] PE (i,2j+1) =cos(i / 10000 2j / e )

[0020] where i∈[1,n] represents the i-th residue, is the dimension, e=128 is the embedding dimension. Add the solved positional encoding to the initial representation E p In the above example, we obtain the position-aware representation:

[0021] S p =E p +PE p

[0022] Then use the corresponding weight transformation matrix W Q , W K and W V is the weight matrix trained by back propagation of the neural network, representing the position perception S p Transformed to That is Q=XW Q ,K=XW K ,V=XW V . Where Q is the query matrix, K is the key matrix, and V is the value matrix. Subsequently, the obtained attention protein representation can be calculated by the following formula:

[0023] S p,attn =MHSA(S p ,Q,K,V)

[0024] MHSA(·) stands for multi-head self-attention, which enhances the ability to capture global features by stacking multiple self-attention blocks in the global aggregation module. The core idea of ​​the self-attention mechanism is that each element calculates the attention weight of all other elements, and then aggregates features based on these weights, so each element aggregates the features of all elements.

[0025] For the drug branch, the same steps as above are performed and the drug sequence Seq obtained in step 1 is d Input embedding layer to generate initial representation Where m is the length of the drug SMILES sequence and e is the embedding size. Calculate the positional encoding PE d and add it to E d , and the location-aware drug feature S d Feed into a stacked self-attention block to generate an attentive drug representation:

[0026] S d,attn =MHSA(S d ,Q,K,V)

[0027] Step 3: Extract residue features within protein pockets

[0028] Because drugs, as small molecules, have a relatively small size compared to the target protein, they usually bind to a local area on the protein surface called a pocket. Therefore, residue features are extracted only when the residue is located in the pocket, and the attention protein feature S obtained in step 1 is p,attn Input to the protein pocket feature extraction module to obtain the protein pocket residue representation S pocket , which can be expressed as:

[0029] S pocket =extract(S p,attn,idx pocket )

[0030] where idx pocket is a set of indices indicating the residues in the pocket. extract() represents the extraction operation, which is based on the residue index idx pocket , from the protein feature S p,attn Note that although only pocket residues are retained, since the self-attention module in step 2 realizes the information aggregation operation between any pair of residues, each pocket residue contains global protein information, including the features of remote residues. The extracted features are used as the initial feature matrix X pocket , that is, X pocket =S pocket .

[0031] Step 4: Cross-modal feature fusion

[0032] Chemical structure and sequence are two different modalities in describing biomolecules. Each modality has its own unique contribution in feature expression, but insufficient interaction between modalities may limit the representation ability. Therefore, a cross-modal feature fusion module is introduced to improve the representation ability of the native single-modal feature fusion module. This is achieved through the multi-head cross attention (MHCA) mechanism.

[0033] a) Cross-modal fusion of target features

[0034] First, the adjacency matrix of the protein pocket obtained in step 1 and the initial feature matrix X obtained in step 3 are combined. pocket The structural feature G is generated by an encoder consisting of three graph convolutional layers and a graph pooling layer. pocket The graph convolution layer and graph pooling layer can be directly called by PyTorch's graph neural network library torch_geometric. Each graph convolution layer is calculated as follows:

[0035]

[0036] Among them A pocket =A pocket +I n , I n represents an identity matrix, where n is the dimension of the matrix. is the adjacency matrix of the self-connected pocket graph, It is its element The diagonal matrix of Represents an element in the normalized adjacency matrix. is the hidden feature matrix of the lth GCN layer, is the activation function, W (l)is a learnable weight matrix. Next, the hidden features go through a multi-layer perceptron (MLP) block and are finally pooled to generate the global pocket map features:

[0037]

[0038] The multilayer perceptron (MLP) is a feedforward artificial neural network that maps a set of input vectors to a set of output vectors. Its concept was created after the backpropagation algorithm was proposed. MLP(·) includes a linear layer, a layer normalization layer with dropout, and an activation function. SumPool(·) represents global sum pooling.

[0039] Then the order feature S obtained in step 3 is pocket As the source of the query matrix Q, and using the graph features G pocket Construct the other two matrices K and V. Then apply multi-head cross attention to obtain enhanced (i.e., sequence modality-involved) graph representation:

[0040] G pocket,attn =MHCA(S pocket ,G pocket ,Q,K,V)

[0041] G pocket,attn It is the final representation of the pocket in the cross-modal fusion stage.

[0042] b) Cross-modal fusion of drug features

[0043] Similar to protein pocket representation, the two modalities of drug molecules are also fused through a multi-head cross-attention mechanism. For drugs, the two modalities are sequential features S d,attn , S d,attn is the self-attention feature generated by the global aggregation module from the drug SMILES obtained in step 2, and the molecular graph Gra d = {V d ,E d}. It is encoded and pooled using the same encoder as the protein:

[0044]

[0045] Among them A d =A d +I n , is the adjacency matrix of the self-connected molecular graph, I n represents an identity matrix, where n is the dimension of the matrix; Representative A d The normalized matrix of It is its element The diagonal matrix of represents the elements in the normalized adjacency matrix, is the hidden feature matrix of the l+1th GCN layer, is the hidden feature matrix of the lth GCN layer, and W (l) is a learnable weight matrix. Next, use the sequential feature S d,attn Construct the query matrix Q, using the graph features G d Construct the other two matrices K and V. Then apply multi-head cross attention to enhance the graph representation:

[0046] G d,attn =MHCA(S d,attn ,G d ,Q,K,V)

[0047] G d,attn It is the final representation of the drug in the cross-modal fusion stage.

[0048] Step 5: Interaction of cross-entity features

[0049] The processing in step 4 enhances the correlation between different modes of each entity (drug / protein). However, there is still a lack of interaction between entities. Representing drugs or proteins in an isolated manner may hinder the model from identifying key atoms of drug molecules or key residues of proteins, thereby damaging the overall performance of the model. Therefore, the present invention further introduces a cross-entity feature interaction module to improve the representation ability of the isolated method. Technically speaking, cross-entity interaction is also achieved through a multi-head cross-attention mechanism. The present invention first uses the pocket feature G obtained in step 4 pocket,attn As the source of Q, use the drug feature G obtained in step 4 d,attn as a source of K and V. Thus, cross-attention to enhanced drug characterization of protein pockets (denoted as R d ) is calculated as follows:

[0050] R d =MHCA(G pocket,attn ,G d,attn ,Q,K,V)

[0051] Conversely, Q can also be derived from drug characteristics G d,attn , K and V can also be derived from G pocket,attn , the enhanced protein pocket representation of cross-attention paired drugs, denoted as R pocket , calculate R pocket The formula is as follows:

[0052] R pocket =MHCA(G d,attn ,G pocket,attn ,Q,K,V)

[0053] Step 6: Multi-level protein feature extraction

[0054] After the above steps, the characteristics of the drug and protein pocket have been generated. p The input is to the multi-level protein sequence encoder introduced in the present invention for further extracting the target protein features. The encoder consists of multiple basic CNN layers and a hierarchical pooling readout module, and a hierarchical strategy is used to integrate outputs of different scales. The CNN operation and output hidden features are expressed as follows:

[0055]

[0056] Where l = 1, 2, 3 represents the lth layer, is the hidden representation corresponding to this layer, is the embedded protein sequence, where the sequence length is truncated or padded to L because CNN requires the input to have a fixed length and L is a hyperparameter that depends on the data sample.

[0057] Next, the hidden representations generated by each CNN layer are read and spliced ​​through the hierarchical pooling readout module to form a multi-level representation of the target protein. The overall formula is as follows:

[0058]

[0059] Where || represents a concatenation operation, W is a learnable weight, and MaxPool(·) is a maximum pooling layer that slides a window (called a pooling window or pooling kernel) on the input feature map and takes the maximum value in each window as the output. This process traverses the entire feature map until it reaches the edge.

[0060] Finally, R d and R generated in step 5 pocket Compositional drug-target complex characterization (which can be considered as local features) and hierarchical protein characterization R p (which can be regarded as global features) are connected for the final prediction, and the formula is as follows:

[0061]

[0062] in Represents the drug-target affinity value finally predicted by the model, R d represents the enhanced drug representation generated in step 5, R pocket represents the enhanced protein pocket representation generated in step 5, R prepresents the hierarchical protein representation generated in step 6. Note that MLP(·) represents three stacked MLP layers. The model is trained by minimizing the mean squared error loss function, as shown below:

[0063]

[0064] Where L mse represents minimizing the mean square error, represents the drug-target affinity value of the nth sample in the dataset predicted by the model, y n is the basic fact of the nth sample in the dataset, and N represents the total number of samples.

[0065] Compared with the prior art, the present invention has the following advantages:

[0066] The present invention uses the self-attention mechanism to enable protein pocket residues to learn global residue features, thereby achieving the learning of global protein features. In addition, the cross-modal feature fusion mechanism is used to enhance the sequence features of protein and drug structure features. While enhancing the cross-modal feature fusion within the entity, the present invention also uses the cross-entity interaction module to identify the key atoms of drug molecules or the key residues of proteins, thereby improving the overall performance of the model. Finally, the present invention further introduces multi-level protein feature extraction to further extract protein features. BRIEF DESCRIPTION OF THE DRAWINGS

[0067] The present invention will be further described below in conjunction with the accompanying drawings.

[0068] Figure 1 is a flow chart of the present invention;

[0069] Figure 2 The PDB IDs of the drug-target complex of the present invention are: 1OPI, 3Q7Q, 2BFQ, and PDB ID verification analysis. DETAILED DESCRIPTION

[0070] The present invention is described in detail below in conjunction with the accompanying drawings.

[0071] Example 1

[0072] like Figure 1 As shown, a drug-target affinity prediction method based on a multi-scale hybrid attention network proposed in the present invention comprises the following steps:

[0073] 1. Construction of drug SMILES and molecular graphs as well as residue sequences and pocket graphs of target proteins

[0074] Drug molecules are usually represented by their natural sequence descriptors SMILES, which can be defined as Seq d ={c1 ,c 2 ,...,c m}, where each element c i is the character that represents the chemical bond or atom in the SMILES string. i There are 63 characters in total, so Seq d It is encoded as a numerical sequence consisting of numbers in [0,63). After obtaining SMILES, the molecular graph is calculated through SMILES. First, the drug SMILES is converted into a molecular graph using RDKit. RDKit is an open source toolkit for chemical informatics. It is based on the operation of 2D and 3D molecules of compounds and uses machine learning methods to generate compound descriptors, fingerprints, calculate compound structural similarity, 2D and 3D molecular displays, etc. It can be called and used based on the Python language; the drug graph descriptor Gra d = {V d ,E d},V d represents the set of nodes containing atoms in a molecule, E d is the bond edge set. In order to generate the initial graph features, each atom is described by a 63-dimensional vector, which represents the chemical properties of the atom. The 63-dimensional vectors are: atom type is represented by a 44-dimensional vector, degree (number of covalent bonds) is represented by a 6-dimensional vector, implicit valence (the implicit valence of the atom) is represented by a 6-dimensional vector, implicit quantity (the number of implicit hydrogen atoms around the atom) is represented by a 6-dimensional vector, aromaticity (whether the atom belongs to the aromatic system) is represented by a 1-dimensional vector, and the total is represented by a 63-dimensional vector. Therefore, the initial feature matrix can be obtained and the adjacency matrix Where N d =|V d | represents the number of atomic nodes in the graph, D d =63 represents the characteristic dimension of the atom, Represents the real number field, indicating that all elements of this matrix are real numbers.

[0075] Proteins are naturally composed of sequences of amino acid residues. Therefore, the target protein can be directly represented as a residue sequence Seq p = {r 1 ,r 2 ,...,r n}, where each element r i Represents a residue. There are 20 types of residues in total, plus one filler value type, so Seq pis encoded as a sequence of numbers in [0,21). For graph construction, in order to better focus on the local binding mode, we construct a graph for the local protein pocket rather than the entire protein. Residues whose distance to the nearest drug molecule atom is less than a certain threshold δ are considered pocket residues. The threshold δ hyperparameter is set to 5 angstroms in the present invention and can be further adjusted to 3 to 10 angstroms based on experimental measurements. Then, the edge set E is generated based on the Euclidean distance between the Cα atom pairs of any two residues. pocket Only when the distance is less than a certain threshold value θ, two residual nodes will be connected. The threshold value θ is a hyperparameter, which is set to 6 angstroms in the present invention and can be further adjusted to 3 to 10 angstroms according to experimental measurements. Based on the above process, the descriptor Gra of the protein pocket can be generated. pocket = {V pocket ,E pocket},V d represents the node set containing the residues in the protein pocket, E d Represents the edge set of the protein pocket residue contact graph. And the adjacency matrix (where N pocket =|V pocket | indicates the number of pocket residues).

[0076] 2. Using multi-head self-attention mechanism to achieve global feature aggregation

[0077] Multi-head self-attention (MHSA) is used in the present invention to process protein sequences and drug SMILES strings. Unlike convolutional networks, the self-attention mechanism enables an element (e.g., a residue in a protein sequence) to perceive information from distant positions, thereby helping to aggregate global representations. This feature is particularly important in the context of the present invention, such as aggregating global features of residues within a binding pocket. Formally, MHSA takes an initial feature X = R n×d , the scaled dot-product attention is calculated as follows:

[0078]

[0079] in d k represents the feature dimension of the query matrix and the key matrix, d v represents the characteristic dimension of the value matrix, K T Represents the transpose of the K matrix. It is the linear transformation of the input residual feature matrix X, that is, Q = XW Q ,K=XW K ,V=XW V , three weights Based on the scaled dot product, multiple heads can be introduced in parallel to learn different attention weights, thereby enriching the learned representation:

[0080] head i =Attention(Q i ,K i ,V i )

[0081] H attn =concat(head i ,head 2 ,…,head h )W O

[0082] in H represents the i-th self-attention head. attn Represents the overall multi-head attention matrix formed by concatenating and mapping h attention heads. concat represents the concatenation operation on the feature dimension.

[0083] The dimensions of the connected h heads are pass Transform it back to the original feature dimension. Then perform residual connection and layer normalization:

[0084] H attn ′=LayerNorm(X+H attn )

[0085] H attn ' represents the overall multi-head attention matrix after residual connection and layer normalization. LayerNorm(·) represents layer normalization, which can be expressed by the following formula where x i is the activation value of the i-th layer, μ and σ 2 are the mean and variance of all activation values ​​of the layer, ε is a very small number to avoid division by zero, and γ and β are learnable parameters. The obtained features are input into the feedforward network (FFN), and then a round of residual connection and layer normalization is performed to generate the final features. The formula is as follows:

[0086] X attn =LayerNorm(H attn ′+FFN(H attn ′)

[0087] To avoid repeating the above process, unless otherwise specified, the present invention expresses MHSA as the following equation:

[0088] X attn =MHSA(X,Q,K,V)

[0089] For the protein branch, the residue sequence Seq obtained in step 1 is p Input embedding layer to generate initial representation Where n is the length of the protein sequence and e = 128 is the embedding dimension. Next, the generated data is encoded with the same dimension through the sinusoidal encoding scheme, i.e. Here is the encoding process:

[0090] PE (i,2j) = sin(i / 10000 2j / e )

[0091] PE (i,2j+1) =cos(i / 10000 2j / e )

[0092] where i∈[1,n] represents the i-th residue, is the dimension, e=128 is the embedding dimension. Add the solved positional encoding to the initial representation E p In the above example, we obtain the position-aware representation:

[0093] S p =E p +PE p

[0094] Then use the corresponding weight transformation matrix W Q , W K and W V is the weight matrix trained by back propagation of the neural network, representing the position perception S p Transformed to That is Q=XW Q ,K=XW K ,V=XW V . Where Q is the query matrix, K is the key matrix, and V is the value matrix. Subsequently, the obtained attention protein representation can be calculated by the following formula:

[0095] S p,attn =MHSA(S p ,Q,K,V)

[0096] MHSA(·) stands for multi-head self-attention, which enhances the ability to capture global features by stacking multiple self-attention blocks in the global aggregation module. The core idea of ​​the self-attention mechanism is that each element calculates the attention weight of all other elements, and then aggregates features based on these weights, so each element aggregates the features of all elements.

[0097] For the drug branch, the same steps as above are performed and the drug sequence Seq obtained in step 1 is d Input embedding layer to generate initial representation Where m is the length of the drug SMILES sequence and e is the embedding size. Calculate the positional encoding PE dand add it to E d , and the location-aware drug feature S d Feed into a stacked self-attention block to generate an attentive drug representation:

[0098] S d,attn =MHSA(s d ,Q,K,V)

[0099] 3. Extracting residue features within protein pockets

[0100] Because drugs, as small molecules, have a relatively small size compared to the target protein, they usually bind to a local area on the protein surface called a pocket. Therefore, residue features are extracted only when the residue is located in the pocket, and the attention protein feature S obtained in step 1 is p,attn Input to the protein pocket feature extraction module to obtain the protein pocket residue representation S pocket , which can be expressed as:

[0101] S pocket =extract(S p,attn ,idx pocket )

[0102] where idx pocket is a set of indices indicating the residues in the pocket. extract() represents the extraction operation, which is based on the residue index idx pocket , from the protein feature S p,attn Note that although only pocket residues are retained, since the self-attention module in step 2 realizes the information aggregation operation between any pair of residues, each pocket residue contains global protein information, including the features of remote residues. The extracted features are used as the initial feature matrix X pocket , that is, X pocket =S pocket .

[0103] 4. Cross-modal feature fusion

[0104] Chemical structure and sequence are two different modalities in describing biomolecules. Each modality has its own unique contribution in feature expression, but insufficient interaction between modalities may limit the representation ability. Therefore, a cross-modal feature fusion module is introduced to improve the representation ability of the native single-modal feature fusion module. This is achieved through the multi-head cross attention (MHCA) mechanism.

[0105] The multi-headed cross-attention (MHCA) mechanism focuses on the interaction between the two sides. This is particularly useful for cross-modal fusion (such as sequential and structural modes) and cross-entity interactions (such as between drugs and proteins). For two modalities (or entities) X and Y, if the query matrix Q is transformed from X and K, and V is transformed from Y, then the cross-attention can be regarded as an enhancement operation on Y, calculated as:

[0106] Y attn =MHCA(X,Y,Q,K,V)

[0107] Where MHCA(·) is the same as MHSA(·) above, except that Q and K, V have different transformation sources. Conversely, if Y is regarded as the source of the query matrix, MHCA(·) can also be used to enhance X:

[0108] X attn =MHCA(X,Y,Q,K,V)

[0109] a) Cross-modal fusion of target features

[0110] First, the adjacency matrix of the protein pocket obtained in step 1 and the initial feature matrix X obtained in step 3 are combined. pocket The structural feature G is generated by an encoder consisting of three graph convolutional layers and a graph pooling layer. pocket The graph convolution layer and graph pooling layer can be directly called by PyTorch's graph neural network library torch_geometric. Each graph convolution layer is calculated as follows:

[0111]

[0112] Among them A pocket =A pocket +I n , I n represents an identity matrix, where n is the dimension of the matrix. is the adjacency matrix of the self-connected pocket graph, It is its element The diagonal matrix of Represents an element in the normalized adjacency matrix. is the hidden feature matrix of the lth GCN layer, is the activation function, W (l) is a learnable weight matrix. Next, the hidden features go through a multi-layer perceptron (MLP) block and are finally pooled to generate the global pocket map features:

[0113]

[0114] The multilayer perceptron (MLP) is a feedforward artificial neural network that maps a set of input vectors to a set of output vectors. Its concept was created after the backpropagation algorithm was proposed. MLP(·) includes a linear layer, a layer normalization layer with dropout, and an activation function. SumPool(·) represents global sum pooling.

[0115] Then the order feature S obtained in step 3 is pocket As the source of the query matrix Q, and using the graph features G pocket Construct the other two matrices K and V. Then apply multi-head cross attention to obtain enhanced (i.e., sequence modality-involved) graph representation:

[0116] G pocket,attn =MHCA(S pocket ,G pocket ,Q,K,V)

[0117] G pocket,attn It is the final representation of the pocket in the cross-modal fusion stage.

[0118] b) Cross-modal fusion of drug features

[0119] Similar to protein pocket representation, the two modalities of drug molecules are also fused through a multi-head cross-attention mechanism. For drugs, the two modalities are sequential features S d,attn , S d,attn is the self-attention feature generated by the global aggregation module from the drug SMILES obtained in step 2, and the molecular graph Gra d = {V d ,E d}. It is encoded and pooled using the same encoder as the protein:

[0120]

[0121]

[0122] Among them A d =A d +I n , is the adjacency matrix of the self-connected molecular graph, I n represents an identity matrix, where n is the dimension of the matrix. It is its element The diagonal matrix of Represents an element in the normalized adjacency matrix. is the hidden feature matrix of the l+1th GCN layer, is the hidden feature matrix of the lth GCN layer, and W (l)is a learnable weight matrix. Next, use the sequential feature S d,attn Construct the query matrix Q, using the graph features G d Construct the other two matrices K and V. Then apply multi-head cross attention to enhance the graph representation:

[0123] G d,attn =MHCA(S d,attn ,G d ,Q,K,V)

[0124] G d,attn It is the final representation of the drug in the cross-modal fusion stage.

[0125] 5. Cross-Entity Feature Interaction

[0126] The processing in step 4 enhances the correlation between different modes of each entity (drug / protein). However, there is still a lack of interaction between entities. Representing drugs or proteins in an isolated manner may hinder the model from identifying key atoms of drug molecules or key residues of proteins, thereby damaging the overall performance of the model. Therefore, the present invention further introduces a cross-entity feature interaction module to improve the representation ability of the isolated method. Technically speaking, cross-entity interaction is also achieved through a multi-head cross-attention mechanism. The present invention first uses the pocket feature G obtained in step 4 pocket,attn As the source of Q, use the drug feature G obtained in step 4 d,attn as a source of K and V. Thus, cross-attention to enhanced drug characterization of protein pockets (denoted as R d ) is calculated as follows:

[0127] R d =MHCA(G pocket,attn ,G d,attn ,Q,K,V)

[0128] Conversely, Q can also be derived from drug characteristics G d,attn , K and V can also be derived from G pocket,attn , the enhanced protein pocket representation of cross-attention paired drugs, denoted as R pocket , calculate R pocket The formula is as follows:

[0129] R pocket =MHCA(G d,attn ,G pocket,attn ,Q,K,V)

[0130] 6. Multi-level protein feature extraction

[0131] After the above steps, the characteristics of the drug and protein pocket have been generated. pThe input is to the multi-level protein sequence encoder introduced in the present invention for further extracting the target protein features. The encoder consists of multiple basic CNN layers and a hierarchical pooling readout module, and a hierarchical strategy is used to integrate outputs of different scales. The CNN operation and output hidden features are expressed as follows:

[0132]

[0133] Where l = 1, 2, 3 represents the lth layer, is the hidden representation corresponding to this layer, is the embedded protein sequence. Here the sequence length is truncated or padded to L, because CNN requires the input to have a fixed length, and L is a hyperparameter that depends on the data sample.

[0134] Next, the hidden representations generated by each CNN layer are read and spliced ​​through the hierarchical pooling readout module to form a multi-level representation of the target protein. The overall formula is as follows:

[0135]

[0136] Where || represents a concatenation operation, W is a learnable weight, and MaxPool(·) is a maximum pooling layer that slides a window (called a pooling window or pooling kernel) on the input feature map and takes the maximum value in each window as the output. This process traverses the entire feature map until it reaches the edge.

[0137] Finally, R d and R generated in step 5 pocket Compositional drug-target complex characterization (which can be considered as local features) and hierarchical protein characterization R p (which can be regarded as global features) are connected for the final prediction, and the formula is as follows:

[0138]

[0139] in Represents the drug-target affinity value finally predicted by the model, R d represents the enhanced drug representation generated in step 5, R pocket represents the enhanced protein pocket representation generated in step 5, R p represents the hierarchical protein representation generated in step 6. Note that MLP(·) represents three stacked MLP layers. The model is trained by minimizing the mean squared error loss function, as shown below:

[0140]

[0141] Where L mse represents minimizing the mean square error, represents the drug-target affinity value of the nth sample in the dataset predicted by the model, y n is the basic fact of the nth sample in the dataset, and N represents the total number of samples.

[0142] Example 2

[0143] The MHAN-DTA method proposed in the present invention has achieved leading performance in multiple benchmark datasets. In this embodiment, the PDBbind-v2020 dataset is divided into training set, validation set, and test set (i.e., internal test set) according to the ratio of 8:1:1. After training the training set and adjusting the hyperparameters of the validation set, the internal test set (inner test) and the external benchmark test set (CASF-2013, CASF-2016, CSAR-HiQ) have achieved international leading levels in the same period. Here, a case analysis on the attention mechanism is used to demonstrate the interpretability of MHAN-DTA. Specifically, in Figure 2 Three representative drug-target complexes are visualized in Figure 1 (PDB ID: 1OPI, 3Q7Q, 2BFQ, PDB ID is a unique identifier for each entry in the Protein Data Bank). The purple residues annotated in the second column are pocket residues with high attention scores and are considered to be very important in forming drug-target complexes. The third column shows the distance of the highest-scoring residues from the center of the drug molecule. Usually the highest-scoring residues are close to the center of the drug molecule and highly overlap with the pocket residues. These observations demonstrate the interpretability of MHAN-DTA.

[0144] Example 3 Discussion of Limitations

[0145] This algorithm is generally applicable to drug-target affinity prediction tasks, and its overall performance has reached the international leading level of the same period. In order to further explore the potential defects of MHAN-DTA, this embodiment collects multiple samples with the worst performance in the internal test set, trying to find out whether they have common characteristics. However, these samples are no different from the best performing group in terms of average protein length, family and function. This shows that there is no significant commonality in the poorly performing samples, and specific analysis is required for specific situations. Here, the present invention lists two representative groups.

[0146] (1) Vascular endothelial growth factor receptor II (VEGFR-2, UniProtIDP35968, UniProtID is the number used to uniquely identify each protein entry in the UniProt database). Among the worst performing samples, there are three samples related to VEGFR-2, namely 3CPB (rmse = 3.409), 3CPC (rmse = 2.647), and 3VHK (rmse = 3.168). rmse is the root mean square error, which is a common indicator to measure the difference between the model prediction value and the true value. It is obtained by averaging the squares of the prediction errors and then taking the square root. The complexity of this protein stems from its complex domain architecture. Specifically, it has up to eight domains, including seven Ig-like C2-type domains and a protein kinase domain. These domains may confuse the model when global features are aggregated, making it unable to capture the most relevant information. In addition, these domains may undergo conformational changes when activated, which may also affect the predictions of MHAN-DTA trained on static protein structures.

[0147] (2) Carbonic anhydrase II (CA2, UniProtIDP00918). Among the worst performing samples, there are four samples related to CA2, namely 4E4A (rmse = 2.775), 5FNJ (rmse = 2.497), 5LLG (rmse = 2.563) and 5VGY (rmse = 2.911). Its active site contains a zinc ion, which is key to its catalytic activity. The highly specialized active site of CA2, including the binding of the zinc ion and its specific orientation with water molecules, constitutes a complex environment that MHAN-DTA may find difficult to replicate. Deep learning models often summarize protein-ligand interactions based on broader features but may fail to capture the nuances of CA2's unique metal coordination environment.

Claims

1. A drug-target affinity prediction method based on a multi-scale hybrid attention network, characterized in that: The following steps are involved: Step 1: Construction of drug SMILES and molecular graph as well as residue sequence and pocket graph of target protein; The drug molecule is represented by its natural sequence descriptor SMILES, which is defined as Seq d ={c1,c2,...,c m }, where each element c i The characters representing chemical bonds or atoms in the SMILES string are used to convert the drug SMILES into a molecular graph and obtain the drug graph descriptor Gram d = {V d ,E d },V d represents the set of nodes containing atoms in a molecule, E d is the key edge set; get the initial feature matrix and the adjacency matrix Where N d =|V d | represents the number of atomic nodes in the graph, D d =63 represents the characteristic dimension of the atom, Represents the real number domain, that is, all elements of this matrix are real numbers; The target protein is represented as a residue sequence Seq p ={r1,r2,...,r n }, where each element r i Represents a residue; constructs a graph for a local protein pocket rather than the entire protein; generates a descriptor Gram for the protein pocket pocket = {V pocket ,E pocket },V d represents the node set containing the residues in the protein pocket, E d represents the contact graph edge set of protein pocket residues, E pocket The Euclidean distance between the Cα atom pairs of any two residues generates an edge set, an adjacency matrix Where N pocket =|V pocket | indicates the number of pocket residues; Step 2: Use multi-head self-attention mechanism to achieve global feature aggregation; The residue sequence Seq obtained in step 1 p Input embedding layer to generate initial representation Where n is the length of the protein sequence and e = 128 is the embedding dimension; Next, the generated data is encoded with the same dimensions using a sinusoidal encoding scheme, i.e. Add the solved position encoding to the initial representation E p , we get the position-aware representation: S p =E p +PE p Then use the corresponding weight transformation matrix Position-aware representations p Transformed to That is Q=XW Q ,K=XW K ,V=XW V ; Q is the query matrix, K is the key matrix, V is the value matrix; W Q , W K and W V is the weight matrix trained by back propagation of the neural network. Then, the obtained attention protein representation is calculated by the following formula: W p,attn =MHSA(S p ,Q,K,V) Where MHSA(·) stands for multi-head self-attention, For the drug branch, the same steps as above are performed and the drug sequence Seq obtained in step 1 is d Input embedding layer to generate initial representation Where m is the length of the drug SMILES sequence and e is the embedding size; calculate the positional encoding PE d and add it to E d , and the location-aware drug feature S d Feed into a stacked self-attention block to generate an attentive drug representation: S d,attn =MHSA(S d ,Q,K,V) Step 3: Extract residue features within protein pockets The attention protein feature S obtained in step 1 p,attn Input to the protein pocket feature extraction module to obtain the protein pocket residue representation S pocket , which can be expressed as: W pocket =extract(S p,attn ,idx pocket ) where idx pocket is a set of indices indicating the residues in the pocket. extract() represents the extraction operation, which is based on the residue index idx pocket , from the protein feature S p,attn The corresponding pocket residue representation is extracted from the matrix X. Although only pocket residues are retained, since the self-attention module in step 2 realizes the information aggregation operation between any pair of residues, each pocket residue contains global protein information, including the features of remote residues. The extracted features are used as the initial feature matrix X. pocket , that is, X pocket =S pocket . Step 4: Cross-modal feature fusion A cross-modal feature fusion module is introduced to improve the representation ability of the native single-modal feature fusion module, which is achieved through the multi-head cross-attention MHCA mechanism, including the cross-modal fusion of target features and the cross-modal fusion of drug features; finally, the final representation G of the pocket in the cross-modal fusion stage is obtained. pocket,attn The final representation of the drug in the cross-modal fusion stage G d,attn ;in, G pocket,attn =MHCA(S pocket ,G pocket ,Q,K,V) S pocket Represent the sequential features as the source of the query matrix Q and use the graph features G pocket Construct the other two matrices K and V; G d,attn =MHCA(S d,attn ,G d ,Q,K,V) S d,attn Represents sequential features to construct the query matrix Q, using graph features G d Construct the other two matrices K and V; Step 5: Interaction of cross-entity features First, use the pocket feature G obtained in step 4 pocket,attn As the source of Q, use the drug feature G obtained in step 4 d,attn Enhanced drug characterization of cross-attention protein pockets as a source of K and V (denoted as R d ) is calculated as follows: R d =MHCA(G pocket,attn ,G d,attn ,Q,K,V) Conversely, Q can also be derived from drug characteristics G d,attn , K and V can also be derived from G pocket,attn , the enhanced protein pocket representation of cross-attention paired drugs, denoted as R pocket , calculate R pocket The formula is as follows: R pocket =MHCA(G d,attn ,G pocket,attn ,Q,K,V) Step 6: Multi-level protein feature extraction The E obtained in step 2 through the embedding layer p The input is sent to the introduced multi-level protein sequence encoder for further extracting the target protein features. The encoder consists of multiple basic CNN layers and a hierarchical pooling readout module, which uses a hierarchical strategy to integrate outputs of different scales. The CNN operation and output hidden features are expressed as follows: Where l = 1, 2, 3 represents the lth layer, is the hidden representation corresponding to this layer, is the embedded protein sequence. Here the sequence length is truncated or padded to L, because CNN requires the input to have a fixed length, and L is a hyperparameter that depends on the data sample; Next, the hidden representations generated by each CNN layer are read and spliced ​​through the hierarchical pooling readout module to form a multi-level representation of the target protein; the overall formula is as follows: Where || represents a concatenation operation, W is a learnable weight, and MaxPool(·) is a maximum pooling layer that slides a window (called a pooling window or pooling kernel) on the input feature map and takes the maximum value in each window as the output. This process traverses the entire feature map until it reaches the edge; Finally, R d and R generated in step 5 pocket Compositional drug-target complex characterization and hierarchical protein characterization p Connect them together for the final prediction, the formula is as follows: in Represents the drug-target affinity value finally predicted by the model, R d represents the enhanced drug representation generated in step 5, R pocket represents the enhanced protein pocket representation generated in step 5, R p represents the hierarchical protein representation generated in step 6. Note that MLP(·) represents three stacked MLP layers. The model is trained by minimizing the mean squared error loss function, as shown below: Where L mse represents minimizing the mean square error, represents the drug-target affinity value of the nth sample in the dataset predicted by the model, y n is the basic fact of the nth sample in the dataset, and N represents the total number of samples.

2. According to claim 1, a drug-target affinity prediction method based on a multi-scale hybrid attention network is characterized in that: In step 1, c i There are 63 characters in total, so Seq d is encoded as a sequence of numbers in the range [0,63); each element r i Represents a residue. There are 20 types of residues in total, plus one filler value type, Seq p is encoded into a sequence of numbers in [0,21); the threshold δ hyperparameter is between 3 and 10 angstroms.

3. According to claim 1, a drug-target affinity prediction method based on a multi-scale hybrid attention network is characterized in that: In step 2, the encoder encoding process: ON (i,2j) =sin(i / 10000 2j / e ) ON (i,2j+1) =cos(i / 10000 2j / e ) where i∈[1,n] represents the i-th residue, is the dimension, e=128 is the embedding dimension; Add the solved position encoding to the initial representation E p In the above example, we obtain the position-aware representation S p .

4. According to claim 1, a drug-target affinity prediction method based on a multi-scale hybrid attention network is characterized in that: In step 2, for the global features of the residues in the polymer binding pocket, formally, MHSA takes an initial feature X = R n×d , the scaled dot-product attention is calculated as follows: in d k represents the feature dimension of the query matrix and the key matrix, d v represents the characteristic dimension of the value matrix, K T Represents the transpose of the K matrix, which is the linear transformation of the input residual feature matrix X, that is, Q = XW Q ,K=XW K ,V=XW V , three weights Based on the scaled dot product, multiple heads are introduced in parallel to learn different attention weights to enrich the learned representation: head i =Attention(Q i ,K i ,V i ) H attn =concat(head i ,head2,…,head h )W O in represents the i-th self-attention head, H attn It represents the overall multi-head attention matrix formed by concatenating and mapping h attention heads, and concat represents the concatenation operation on the feature dimension. The dimensions of the connected h heads are pass Transform it back to the original feature dimension. Then perform residual connection and layer normalization: H attn ′=LayerNorm(X+H attn ) H attn ' represents the overall multi-head attention matrix after residual connection and layer normalization, LayerNorm(·) represents layer normalization, which is expressed by the following formula where x i is the activation value of the i-th layer, μ and σ 2 are the mean and variance of all activation values ​​of the layer, ε is a very small number to avoid division by zero, γ and β are learnable parameters; the obtained features are input into the feedforward network FFN, and then a round of residual connection and layer normalization are performed to generate the final features. The formula is as follows: X attn =LayerNorm(H attn ′+FFN(H attn ′) Express MHSA as the following equation: X attn =MHSA(X,Q,K,V) 5. The drug-target affinity prediction method based on a multi-scale hybrid attention network according to claim 1, characterized in that: In step 4, the cross-modal fusion steps of target features are as follows: First, the adjacency matrix of the protein pocket obtained in step 1 and the initial feature matrix X obtained in step 3 are combined. pocket The structural feature G is generated by an encoder consisting of three graph convolutional layers and a graph pooling layer. pocket , where the graph convolution layer and graph pooling layer are directly called by PyTorch's graph neural network library torch_geometric, and each graph convolution layer is calculated as follows: Among them A pocket =A pocket +I n , I n represents an identity matrix, where n is the dimension of the matrix. is the adjacency matrix of the self-connected pocket graph, It is its element The diagonal matrix of Represents the elements in the normalized adjacency matrix; is the hidden feature matrix of the lth GCN layer, is the activation function, W (l) is a learnable weight matrix. Next, the hidden features go through a multi-layer perceptron (MLP) block and are finally pooled to generate the global pocket map features: MLP(·) includes a linear layer, a layer normalization layer with dropout, and an activation function. SumPool(·) represents global sum pooling. Then the order feature S obtained in step 3 is pocket As the source of the query matrix Q, and using the graph features G pocket Construct the other two matrices K and V, and then apply multi-head cross attention to obtain the enhanced graph G pocket,attn .

6. The drug-target affinity prediction method based on a multi-scale hybrid attention network according to claim 1, characterized in that: In step 4, the cross-modal fusion of drug features is as follows: For drugs, the two modes are the sequential features S d,attn , S d,attn is the self-attention feature generated by the global aggregation module from the drug SMILES obtained in step 2, and the molecular graph Gra d = {V d ,E d }; encoded and pooled using the same encoder as the protein: Among them A d =A d +I n , is the adjacency matrix of the self-connected molecular graph, I n represents an identity matrix, where n is the dimension of the matrix; It is its element The diagonal matrix of Represents the elements in the normalized adjacency matrix; is the hidden feature matrix of the l+1th GCN layer, is the hidden feature matrix of the lth GCN layer, and W (l) is a learnable weight matrix; next, use the sequential feature S d,attn Construct the query matrix Q, using the graph features G d Construct the other two matrices K and V. Then apply multi-head cross-attention enhanced graph representation G d,attn .

Images