Virtual screening method based on topological pharmacophore molecular graph modeling and target feature analysis

Abstract

The present application belongs to the technical field of virtual screening, and particularly relates to a virtual screening method based on topological pharmacophore graph modeling and target feature analysis. The present application performs virtual screening through topological pharmacophore graph modeling and target position score matrix analysis to predict the interaction between a given drug and a target. First, key pharmacophore atoms in a molecule are found, a topological pharmacophore graph is established, an algorithm for graph reduction is designed to obtain a reduced graph without changing the original topological structure, and then graph neural network is used to learn the topological knowledge and pharmacophore knowledge to obtain the feature representation of the drug molecule. Second, the position score matrix of the target is analyzed from the angles of global, adjacent and local fragments to learn the potential feature information of the amino acid sequence to obtain the feature representation of the target protein. Then, the feature representations of the drug and the target are learned based on the method of multi-kernel learning to obtain the interaction feature. Finally, the interaction between the drug and the target is predicted based on a multi-layer perception machine.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of virtual screening, and particularly relates to a virtual screening method based on topological pharmacophore molecular graph modeling and target feature analysis. BACKGROUND

[0002] Virtual screening (VS) is an important part of the new drug development process, which uses computer technology to search and filter candidate compounds from large compound databases. VS greatly improves the speed of new drug development and reduces costs compared to traditional experimental research methods.

[0003] With the development of artificial intelligence technology and the increase of available chemical and biological activity data, AI has been widely used in VS. VS can be roughly divided into two categories: ligand-based virtual screening and structure-based virtual screening. One of the typical directions of ligand-based virtual screening methods is to predict the interaction between drugs and targets. The specific implementation is mainly based on two ideas: first, analyzing similarity information, that is, starting from the similarity information of drugs and the similarity information of targets, capturing the relationship of interaction from the similarity information through matrix decomposition, least squares method, gradient boosting machine, etc., and finally predicting the interaction between drugs and targets. Second, extracting features, that is, respectively representing drugs and target proteins, taking the descriptors of drugs and targets as input, extracting potential feature information, and finally predicting the interaction based on the feature representation of both.

[0004] The descriptors of drug molecules can be divided into one-dimensional, two-dimensional, three-dimensional, etc. according to the required molecular structure dimension for calculation, among which the representative ones are one-dimensional SMILES sequence, two-dimensional molecular substructure fingerprint, and three-dimensional molecular crystal structure. Most targets are protein macromolecules, and their primary sequences (amino acids), motifs, secondary structures, and folds are generally analyzed.

[0005] The above technology lays the foundation for the research of virtual screening, but has the following shortcomings: first, VS based on similarity information mainly starts from known drugs and targets, and the generalization of the model is low when facing new sample data, which is not conducive to the discovery of new drugs, and the method has low interpretability. Second, most of the current VS based on features starts from one-dimensional SMILES sequence and amino acid sequence, and the feature information extracted from one-dimensional sequence is limited; some studies collect the chemical attribute information of each atom by establishing a molecular graph and learning the features of the molecular graph through a graph neural network, which is difficult to reflect the influence of key pharmacophores in a molecular structure. Third, the current background technology does not fully represent the target, and the information is too single. Moreover, the drug and target features are generally directly spliced for feature fusion, which is too simple. SUMMARY

[0006] In view of the above problems, the present application provides a virtual screening method based on topological pharmacophore molecular graph modeling and target position-specific score matrix analysis. The key pharmacophore molecular features are extracted through the topological pharmacophore molecular graph. The position-specific matrix is calculated according to the amino acid sequence of the target protein, and the potential pattern sequence feature information is extracted through various feature analysis. Based on the above two types of feature information, the prediction and screening classification work for new drug research and development is used.

[0007] The technical scheme of the present application is:

[0008] The virtual screening method based on topological pharmacophore molecular graph modeling and target feature analysis comprises:

[0009] The SMILES formula of the drug molecule and the amino acid sequence of the target protein are extracted respectively, specifically:

[0010] Based on the topological pharmacophore molecular graph, the molecular graph features of the topological pharmacophore molecular graph are extracted by using the graph neural network. The construction method of the topological pharmacophore molecular graph is: according to the input SMILES formula of the drug molecule, the pharmacophore atoms in the molecule and their indexes are analyzed by using the RDkit tool, and the adjacency matrix of the molecule is calculated, and then the molecular graph is established;

[0011] Based on the position-specific score matrix analysis, the amino acid sequence of the target protein is extracted. The method is: first, based on the input amino acid sequence of the target protein, the position-specific score matrix is obtained. Each item in the matrix reflects the frequency of each amino acid of the target protein at a specific position. Then, the matrix feature analysis is carried out, and the long-term sequence feature extraction, the adjacent coupling relationship feature extraction and the short-term fragment feature extraction are carried out respectively. The three extracted features are spliced to obtain the amino acid sequence feature of the target protein;

[0012] The molecular graph features and the amino acid sequence features of the target protein are fused to obtain the interaction features;

[0013] The interaction features are input into the multilayer perceptron, and the interaction between the drug and the target is predicted after training.

[0014] The present application has the following advantages:

[0015] Compared with the ordinary two-dimensional molecular graph, the topological pharmacophore molecular graph of the present application not only reflects the topological information, but also reflects the pharmacophore and molecular activity information to a certain extent, which is of great significance for predicting the interaction between the drug and the target, and supplements the field knowledge. The graph reduction algorithm reduces the complexity of the molecule, improves the calculation efficiency and the robustness of the model.

[0016] The position-specific score matrix is analyzed from the perspective of global, adjacent and local fragments to obtain rich feature information, which has important reference significance for downstream tasks.

[0017] The method based on multi-core learning learns the feature representation of drugs and targets, and extracts interaction features, which is beneficial to improve the accuracy of virtual screening prediction. BRIEF DESCRIPTION OF DRAWINGS

[0018] Figure 1 is the overall flowchart of the present application.

[0019] Figure 2 is the model diagram of the design scheme of the present application.

[0020] Figure 3 is the flowchart of the topological pharmacophore graph construction.

[0021] Figure 4 is the visual display of the graph reduction case.

[0022] Figure 5 is the flowchart of the graph reduction algorithm. DETAILED DESCRIPTION

[0023] The present application will be described in detail below with reference to the accompanying drawings.

[0024] The present application takes the SMILES formula of the drug molecule and the amino acid sequence of the target protein as input, and the overall flowchart of the scheme is as shown in Figure 1 .

[0025] 1. Topological pharmacophore graph construction

[0026] The pharmacophore refers to the "pharmacophore characteristic element" that plays an important role in the activity of the drug molecule. These "pharmacophore characteristic elements" are active sites when ligands and receptors interact, which can be specific atoms or atom groups, such as oxygen atoms, hydroxyl groups, carbonyl groups, etc. The topological pharmacophore graph modeling of the present application adopts six common pharmacophoric characteristics, namely hydrogen bond donor (Hydrogen Bond Donor), hydrogen bond acceptor (Hydrogen Bond Acceptor), positive ion (Positively Ionizable), negative ion (Negatively Ionizable), zinc binder (Zinc Binders), and aromatic ring (Aromatic Rings).

[0027] Given a drug molecule, to build a topological pharmacophore molecular graph, it is necessary to find and locate its pharmacophore. The present application combines the RDkit tool to assist in finding the pharmacophore in the molecule. Using the RDkit tool to search for a pharmacophore requires manually importing a feature library, which is essentially a feature definition file (Feature Definition File, FDef) that is manually written according to the pharmacophore to be searched and contains definitions and explanations of atomic types. The atomic type defines the alias of the query pharmacophore atom, and a SMARTS expression defines the atom matched by the feature type, for example, a hydrogen bond acceptor can be defined as follows: “AtomType NAcceptor[n;+0;!X3;!$([n;H1](cc)cc)]”.

[0028] The flowchart of the topological pharmacophore molecular graph construction is shown in Figure 3 The SMILES formula of the drug molecule is taken as the input, and RDkit traverses the atoms according to the rules defined in the feature library and calculates the chemical features, thereby finding the pharmacophore atoms and their indices in the molecule. Combined with the SMILES formula of the molecule, the adjacency matrix of the molecule can be calculated by calling the relevant functions of RDkit, and then the molecular graph can be established. By marking the pharmacophore atoms and their indices on the molecular graph, the topological pharmacophore molecular graph can be obtained. As shown in Figure 2 (a), the pharmacophore of the molecule with id “11712649” is extracted and the topological graph structure is created, and the pharmacophore atoms are marked with corresponding gray scale and identifiers.

[0029] 2. Topological pharmacophore molecular graph reduction

[0030] The molecular graph marked with pharmacophore atoms has certain redundancy, and the original pharmacophore molecular graph is simplified under the premise of ensuring the general topology. The present application designs a graph reduction algorithm and defines two reduction rules.

[0031] (1) Reduction of ordinary atoms

[0032] The purpose of reduction is to find a set of ordinary atoms (non-pharmacophore atoms) that can be combined, which satisfy the following conditions: 1) all are ordinary atoms; 2) the total number of atoms in the set is not less than 2; 3) at most only one atom in the set has a degree greater than 2. The ordinary atom reduction process is shown in Figure 3 The input is the topological pharmacophore molecular graph, and all atoms are traversed one by one. If the currently traversed atom is an ordinary atom, then the neighbor atoms are traversed from this atom, and the traversal process is iterated repeatedly until a set of atoms that satisfy the above merging rules is found. The node features of the merged set are updated, and the atoms and edges that are reduced are deleted, thereby completing the reduction of ordinary atoms. Figure 4 (a) and 4(b) are schematic diagrams of the reduction process of ordinary connecting atoms, Figure 5(a) is the reduction process of common atoms.

[0033] (2) Reduction of only aromatic atom groups

[0034] The purpose of the reduction is to find a group of aromatic atoms that meet the following conditions: 1) all are aromatic atoms; 2) in the same aromatic ring; 3) the atoms are only aromatic atoms and do not contain other pharmacophoric characteristics. As shown in Figure 4 (c), the atoms that meet the same aromatic ring are merged. The aromatic atom reduction process is shown in Figure 5 (b), traverse the atom group of all aromatic rings, if the atoms in the group are only aromatic atoms, the group is considered as a group of atoms that can be merged, the node features of the merged group are updated, and the reduced atoms and edges are deleted, thereby completing the reduction of the aromatic atom group.

[0035] (3) Initial embedding generation of each atom in the molecule

[0036] In constructing the original pharmacophore molecule graph, the initial embedding of each atom in the molecule graph needs to be initialized for subsequent learning and message passing of the molecule graph. The present application generates the initial embedding of the atom based on the hash algorithm, and the atoms with pharmacophoric characteristics are sorted into 32-dimensional feature vectors with uniqueness by the hash algorithm, and the common atoms without pharmacophoric characteristics are also generated in a similar manner. Traverse each atom node in the graph, calculate the initial embedding, if a certain atom node has multiple pharmacophoric characteristics, directly sum to obtain the initial embedding.

[0037] In the graph reduction process, after finding a group of atoms that can be merged, the embeddings of all atoms in the group are summed to obtain a new node embedding. The reduced atoms and edges are deleted to obtain the reduced pharmacophore molecule graph.

[0038] 3. Graph neural network learning of topological pharmacophore molecule graph features

[0039] (1) Representation of topological pharmacophore molecule graph

[0040] Define the molecule graph G=(V, E), V is the set of atoms, V={v t} (t＝1:n) , n represents the number of atoms in the molecule; wherein d v represents the feature dimension of the atom, E is the edge set between atoms, E={(s r ,d r ) (r＝1:m) , m represents the number of chemical bonds in the molecule, s r and d r are the starting atom index and ending atom index of the rth edge, respectively.

[0041] (2) Graph neural network (GNN) learning

[0042] GNN takes node features, edges, and graph structure as input, and learns the representation vector of nodes, where the representation vector of node v is denoted as h v Given node v, the neighbor feature aggregation in the kth iteration is shown in equation (1), the update of the current node feature of node v is shown in equation (2), and finally the feature representation of the entire graph is obtained by aggregating the node features through the readout function as shown in equation (3).

[0043]

[0044]

[0045]

[0046] is the neighborhood of v, AGGREGATE(.) is an aggregation function used to pass and aggregate the information of the neighbor nodes of node v. COMBINE(.) is an update function that updates the neighborhood information passed to the current node. READOUT(.) is a readout function that aggregates all node features to obtain the feature of the entire graph.

[0047] 4. Position-specific score matrix analysis

[0048] (1) Calculate the position-specific score matrix

[0049] The position-specific score matrix (PSSM) is also known as the position weight matrix. Each item in the matrix reflects the frequency of each amino acid of the target protein at a specific position. Effective analysis of PSSM can discover potential feature information of proteins.

[0050] First, based on the ncbi-blast tool, combined with the protein comparison database swissprot, the PSSM is calculated, PSSM ∈ R (L ×20) The matrix form is shown in equation 4, L is the length of the amino acid sequence, 20 reflects the 20 basic amino acids that make up the protein, p i,j represents the score of the amino acid residue in the j-type amino acid replacing the i-type amino acid in the multiple sequence evolution.

[0051]

[0052] (2) Matrix feature analysis

[0053] The sigmoid function is used to map the values of the matrix to the range of 0 to 1, as shown in equation 5. The sigmoid calculation does not change the shape of the matrix.

[0054] PSSM = 1 / (1 + e (-1×PSSM) )#(5)

[0055] (a) long-term sequence feature extraction. Considering the overall features of the sequence, the average value of each column of PSSM can be obtained, and a 20-dimensional vector can be embedded into a protein sequence, as shown in equations 6 and 7:

[0056] X = (x1, x2, x3,.., x 20 ) T #(6)

[0057]

[0058] T represents transposition;

[0059] (b) adjacent coupling feature extraction. Considering the position feature information of adjacent amino acids, the product of the i-th amino acid and the j-th amino acid of the adjacent two rows is summed and averaged, and a 400-dimensional vector can be embedded into a protein sequence, as shown in equations 8 and 9:

[0060] Y = (y 1,1 ,…,y 1,20 ,…,y 2,1 ,…,y 2,20 ,…,y 20,1 ,…,y 20,20 ) T #(8)

[0061]

[0062] (c) short-term fragment feature extraction. Considering the short-term position dependence of the amino acid sequence and the potential feature information of the fragment, five windows are set, each window covers 25% of the next window to the right, and five windows overlap after overlapping to cover the entire length of the amino acid sequence. Each column is averaged in each window, and a 100-dimensional vector is obtained by five windows, as shown in equations 10, 11 and 12:

[0063]

[0064]

[0065] Z = Z1||Z2||Z3||Z4||Z5#(12)

[0066] win nums is the number of windows set, Z q is a 20-dimensional vector obtained by averaging in the q-th window.

[0067] Finally, the feature vector F PSSM= X||Y||Z, is a 520-dimensional vector.

[0068] 5. Multi-kernel learning captures interaction features

[0069] Under the multi-kernel framework, the high-dimensional feature representation of drugs and targets is learned, the basic kernel and weight coefficient are selected, and multiple kernels are combined to capture the potential interaction mode between drugs and targets, and then the interaction features are obtained.

[0070] 6. Multi-layer perception predicts drug and target interaction

[0071] The interaction features are input into the multi-layer perception, the interaction between drugs and targets is predicted by training, and the loss function is shown in formula (13).

[0072]

[0073] Where, y I is the true interaction intensity value, is the predicted interaction intensity value, and N is the number of sample sets.

[0074] The present application predicts the interaction between a given drug and a target by topological pharmacophore molecular graph modeling and target position score matrix analysis. First, the key pharmacophore atoms in the molecule are found, the topological pharmacophore molecular graph is established, the graph reduction algorithm is designed to obtain the reduced graph without changing the original topological structure, and the graph neural network is used to learn the topological knowledge and pharmacophore knowledge to obtain the feature representation of the drug molecule. Second, the position score matrix of the target is analyzed from the angles of the global, adjacent and local fragments to learn the potential feature information of the amino acid sequence and obtain the feature representation of the target protein. Then, the feature representation of the drug and the target is learned based on the multi-kernel learning method to obtain the interaction feature. Finally, the interaction between the drug and the target is predicted based on the multi-layer perception.

Claims

1. A virtual screening method based on topological pharmacophore molecular graph modeling and target feature analysis, characterized in that, Comprise: The SMILES of the drug molecule and the amino acid sequence of the target protein are respectively extracted, specifically: Based on the topological pharmacomolecule graph, the molecular graph features of the topological pharmacomolecule graph are extracted by using the graph neural network; the construction method of the topological pharmacomolecule graph is: according to the input SMILES of the drug molecule, the pharmacophore atoms in the molecule and their indexes are analyzed by using the RDkit tool, and the adjacency matrix of the molecule is calculated, and then the molecular graph is established; Based on the position-specific scoring matrix analysis, the amino acid sequence of the target protein is extracted, the method is: first, based on the input amino acid sequence of the target protein, the position-specific scoring matrix is obtained, each item in the matrix reflects the frequency of each amino acid of the target protein at a specific position; then the matrix feature analysis is carried out, and the long-term sequence feature extraction, the adjacent coupling relationship feature extraction and the short-term fragment feature extraction are carried out respectively, and the extracted three kinds of features are spliced to obtain the amino acid sequence feature of the target protein; The molecular graph features and the amino acid sequence features of the target protein are fused to obtain the interaction features; The interaction features are input into the multilayer perception machine, and the interaction between the drug and the target is predicted after training.

2. The virtual screening method based on topological pharmacophore molecular graph modeling and target feature analysis according to claim 1, characterized in that, After constructing the topological pharmacomolecule graph, the topological pharmacomolecule graph is reduced, specifically: Two kinds of reduction rules are defined as: The first reduction rule: all atoms in the group are ordinary atoms, the total number of atoms in the group is not less than 2, and the maximum number of atoms in the group is only one, and the degree of the atom is greater than 2, then the ordinary atom is reduced to one atom; the definition of the ordinary atom is a non-pharmacophore atom; The second reduction rule: all atoms in the group are aromatic atoms, the atoms in the group are in the same aromatic ring, each atom in the group is only an aromatic atom and does not contain other pharmacophore characteristics, then the atoms in the group are reduced to an aromatic atom; The reduction process is: based on the above two reduction rules, all atoms in the topological pharmacomolecule graph are traversed one by one, if the current traversed atom is an ordinary atom, the neighbor atoms are traversed from the atom, the traversal process is iterated repeatedly, until a set of atoms satisfying the first reduction rule is found, the node features of the combined group are updated, and the reduced atoms and edges are deleted, thereby completing the reduction of the ordinary atoms; if the current traversed atom is an aromatic atom, the atom group of the aromatic ring where the aromatic atom is located is traversed, if all atoms in the group are only aromatic atoms, the group is regarded as a group of atoms that can be combined, the node features of the combined group are updated, and the reduced atoms and edges are deleted, thereby completing the reduction of the aromatic atom group. 3.The virtual screening method based on topological pharmacophore molecular graph modeling and target feature analysis according to claim 2, characterized in that, The method for extracting the molecular graph features of the topological pharmacomolecule graph by using the graph neural network is: Let G = (V, E) be a molecule graph, V is a set of atoms, V = {v t} (t＝1:n) , n represents the number of atoms in the molecule; wherein d v represents the feature dimension of the atom, E is a set of edges between atoms, E = {(s r ,d r )} (r＝1:m) , m represents the number of chemical bonds in the molecule, s r and d r are the starting atom index and the ending atom index of the rth edge, respectively; A graph neural network takes node features, edges, and graph structure as input, and learns a representation vector for each node, where the representation vector for node v is denoted by h v ; given node v, the neighbor feature aggregation at iteration k is given by: The current node feature of the updated node is: Finally, the readout function is used to obtain the features of the whole graph by aggregating the node features: wherein, is the neighborhood of v, AGGREGATE(.) is an aggregation function to pass and aggregate the information of the neighboring nodes of a node v, COMBINE(.) is an update function to update the passed neighborhood information to the current node, and READOUT(.) is a readout function to aggregate all the node features to get the features of the whole graph. 4.The virtual screening method based on topological pharmacophore molecular graph modeling and target feature analysis according to claim 1, wherein, The position specific scoring matrix is calculated as follows: based on the ncbi-blast tool, in combination with the protein comparison database swissprot, a position specific scoring matrix PSSM ∈ R (L×20) : where L is the length of the amino acid sequence, 20 reflects the 20 basic amino acids that make up proteins, p i,j represents the score of an amino acid residue in a multiple sequence evolution in which the amino acid residue of type i is replaced by an amino acid residue of type j.

5. The virtual screening method based on topological pharmacophore molecular graph modeling and target feature analysis according to claim 4, characterized in that, The specific method for matrix feature analysis is: The values of the matrix PSSM are mapped to 0-1 by using the sigmoid function: PSSM = 1 / (1+e (-1×PSSM) ) 1) Long-term sequence feature extraction, average each column of PSSM, and embed a protein sequence into a 20-dimensional vector: X = (x1, x2, x3,.., x 20 ) T T represents transposition; 2) Adjacent coupling feature extraction, the product of the i-th and j-th amino acids of the adjacent two rows is summed and averaged, and a protein sequence is embedded into a 400-dimensional vector: Y = (y 1,1 ,…,y 1,20 ,…,y 2,1 ,…,y 2,20 ,…,y 20,1 ,…,y 20,20 ) T 3) Short-term fragment feature extraction, five windows are set, each window covers 25% of the next window to the right, and the five windows overlap to cover the entire length of the amino acid sequence, the average of each column is calculated in each window, and a 100-dimensional vector is obtained by combining the five windows: Z = Z1||Z2||Z3||Z4||Z5 win nums is the number of windows set, Z q is the 20-dimensional vector obtained by averaging within the qth window; The three extracted features are spliced to obtain a feature vector F PSSM = X||Y||Z, which is a 520-dimensional vector. 6.The virtual screening method based on topological pharmacophore molecular graph modeling and target feature analysis according to claim 5, characterized in that, The method of fusing interaction features is: Under the multi-core framework, learn the high-dimensional feature representation of drugs and targets, select the basic kernel and weight coefficient, and combine multiple kernels to capture the potential interaction patterns between drugs and targets, and then obtain the interaction features.

7. The virtual screening method based on topological pharmacophore molecular graph modeling and target feature analysis according to claim 6, characterized in that, The loss function used to train the multi-layer perception by inputting the interaction features is: where y I is the true interaction strength value, is the predicted interaction strength value, and N is the number of samples in the set.

Images