Protein complex model global quality evaluation method based on deep fusion network

Abstract

A protein complex model global quality evaluation method based on a deep fusion network belongs to the field of bioinformatics and computer application, and comprises the following steps: firstly, setting specific thresholds for maximum sequence redundancy, resolution and residue length in a protein database RCSB PDB to screen out a natural protein structure; for each protein structure, generating a disturbance structure model by using an HDck and xTrimoMultimer docking algorithm and combining a chain orientation disturbance strategy; then, extracting evolutionary features, physical features and geometric features based on each protein structure model, and constructing a deep fusion network fusing a deep convolutional neural network, an axial attention mechanism network and a graph attention network so as to fully capture local and global topological information of a protein compound; and finally, a global quality score is obtained through prediction of a decoding module. According to the method, the global topological information of the protein complex model can be effectively captured, so that the quality of the protein complex model is evaluated.

Description

Technical Field

[0001] The present invention belongs to the fields of bioinformatics and computer applications, and in particular relates to a global quality assessment method for protein complex models based on a deep fusion network. Background Art

[0002] Proteins are crucial molecules in living organisms, involved in virtually all biological processes. They are not only the fundamental building blocks of cells and tissues, but also perform multiple functions, including regulation, catalysis, and information transmission. The function and activity of a protein are directly determined by its three-dimensional structure. Accurately predicting protein structures enables a deeper understanding of their functions, interactions, and mechanisms, thereby advancing progress in various fields. However, the accuracy and reliability of predictions are crucial for subsequent applications. By assessing the accuracy, stability, and physical plausibility of structures, structural quality assessment not only supports protein model identification and refinement, but also supports downstream drug development and disease research. Global protein structure quality assessment is an integral part of protein structure prediction. It helps detect globally unreasonable conformations or folding errors, assesses the stability of the entire protein and its adaptability in vivo, and ensures that the model conforms to known biophysical principles. This not only helps scientists determine the accuracy and reliability of prediction models but also provides strong support for subsequent experimental validation, drug design, disease research, and protein engineering. Summary of the Invention

[0003] In order to overcome the problems of existing global quality assessment of protein structure, the present invention proposes a global quality assessment method for protein complex models based on deep fusion network, and designs a new method for global quality assessment of protein structure. Features are extracted from the input protein structure and sequence as the input of the deep fusion network, capturing the intrinsic characteristics of the protein based on evolution, physics, and geometry. Finally, the global quality score of the protein structure is calculated, aiming to evaluate the rationality, stability and physicochemical feasibility of the overall protein structure.

[0004] The technical solution adopted by the present invention to solve its technical problem is:

[0005] A global quality assessment method for protein complex models based on a deep fusion network. First, specific thresholds are set for maximum sequence redundancy, resolution, and residue length in the protein database (RCSB PDB) to screen out natural protein structures. For each protein structure, the HDock and xTrimoMultimer docking algorithms are used in combination with a chain orientation perturbation strategy to generate a perturbed structural model. Then, evolutionary, physical, and geometric features are extracted based on each protein structure model, and a deep fusion network is constructed that integrates a deep convolutional neural network, an axial attention mechanism network, and a graph attention network to fully capture the local and global topological information of the protein complex. Finally, the global quality score is predicted by a decoding module.

[0006] Furthermore, the method comprises the following steps:

[0007] 1) Extract evolutionary features based on protein models. The process is as follows:

[0008] 1.1) Protein sequence alignment substitution matrix BLOSUM62 matrix;

[0009] 1.2) Protein structure sequence spectrum characteristics;

[0010] 1.3) Calculate the normalized distance between the amino acid and the end of the sequence;

[0011] 2) Extract physical features based on protein models;

[0012] 3) Extracting geometric features based on protein models;

[0013] 4) Producing label data: The label value is the template modeling score TM-score obtained by comparing each structure with the corresponding natural structure. The label value ranges from 0 to 1, where the label value corresponding to the natural structure is 1;

[0014] 5) After concatenating and dimensionally processing the extracted features, they are input into the deep fusion network. The Adam optimizer is used to optimize the network weight parameters and the mean square error (MSE) loss function for training. Finally, the global score of the protein structure is calculated by the decoding module of the deep fusion network to further evaluate the global quality of the protein model.

[0015] Further, the process of 1.1) is:

[0016] 1.1.1) Perform multiple sequence alignment on protein sequences and count the substitution frequencies of each pair of amino acids;

[0017] 1.1.2) Convert the frequency of each pair of amino acid substitutions into probabilities, assuming that the probability of substitution is P(X,Y), where X and Y are two amino acid types;

[0018] 1.1.3) Apply the logarithmic function to transform the substitution probability of each pair of amino acids and calculate their logarithmic ratio:

[0019]

[0020] 1.1.4) Fill the calculated logarithmic ratio values ​​of all amino acid pairs into the matrix to obtain the final BLOSUM62 matrix.

[0021] Furthermore, the process of 1.2) is as follows:

[0022] 1.2.1) Using each protein complex model as a benchmark, use the Foldseek template search tool to search for template structures with a similarity score > 0.5 in the PAcluster80 template library, obtaining a total of n sets of matching templates;

[0023] 1.2.2) Perform structural alignment of the searched template structures using US-align, inserting gap marks ("-") for residue positions that do not form valid atomic correspondences;

[0024] 1.2.3) Based on C β Atomic spatial coordinate information is used to calculate the Euclidean distances of all valid residue pairs between the target structure and each template structure. For residue sites with gaps ("-"), the distance calculation of the position is automatically skipped, and finally a distance distribution matrix containing n+1 groups of data is generated;

[0025] 1.2.4) The distance range is divided into 36 equally spaced intervals (step length ), the probability distribution of the same position in each distance map falling into each interval is counted, and the distance value corresponding to the probability peak is selected as the structural sequence spectrum characteristic parameter of the site.

[0026] Furthermore, the process of 2) is as follows:

[0027] 2.1) One-hot encoding is used to vectorize the amino acid sequence of the protein complex, generating a matrix containing 20 independent binary feature bits, each corresponding to a standard amino acid type;

[0028] 2.2) Extract the main chain torsion angle, and ψ, and calculate the main chain C α -C, NC α , C-N bond length and C α -CN、NC α -C, CNC α bond angle;

[0029] 2.3) Use pyRosetta to extract residue intra-energy features, including intra-residue conformational preferences (p_aa_pp), stereochemical constraints (rama_prepreo, omega), and side chain conformational free energy (fa_dun);

[0030] 2.4) Use pyRosetta to extract the interaction energies between residues, including van der Waals attraction (fa_atr) and repulsion (fa_rep), generalized solvation effect (fa_sol), Lennard-Jones spherical core weighted potential (lk_ball_wtd), charge interaction (fa_elec), main chain-side chain hydrogen bond network (hbond_bb_sc) and side chain hydrogen bond cooperativity (hbond_sc);

[0031] 2.5) Analyze the secondary structure characteristics of protein complexes based on the DSSP algorithm;

[0032] 2.6) Construct a two-dimensional energy distribution matrix based on the distance between residues (C β -C β ) and relative orientation angle as coordinate axes, and the normalized energy distribution heat map is generated by Gaussian kernel density estimation;

[0033] 2.7) Extract the physicochemical properties of amino acids, including isoelectric point, polarity index, acidity / alkalinity, hydrogen bond donor / acceptor capacity, hydrophilicity / hydrophobicity, topological polar surface area, and secondary structure tendency;

[0034] 2.8) Use the Voronota tool to extract the contact area and direction features between residues of the protein model based on Voronoi tessellation.

[0035] Furthermore, the process of 3) is as follows:

[0036] 3.1) Calculate the spatial voxel features of the protein structure: First, define a local coordinate system based on the positions of the main chain atoms with the current residue as the center, where the axis direction is determined by the vector pointing to the N atom of the current residue, the starting point is the C atom of the previous residue, and for the first residue, the axis direction points to the C atom of the current residue. α Atoms; calculate three orthogonal basis vectors by vector product and ensure that the relative position in the local coordinate system is fixed so that it is always in the axis direction and in the axis plane; then, select a fixed coordinate point As the origin, construct a coordinate system with a side length of cube, and further divided into 20×20×20 The number of atoms in each voxel is counted and its center coordinates are calculated. At the same time, a value is assigned to each voxel to characterize the spatial density and distribution characteristics of the area.

[0037] 3.2) Calculate the Euler angle, using C α -C β The vector is the rotation axis, and its value is calculated relative to the reference coordinate system, i.e., NC α -C β α, β, and γ rotation angles in the direction of the plane normal;

[0038] 3.3) Calculate the positional relationship between residues, including: C α -C β -C β -C α The dihedral angle between NC α -C β -C β The dihedral angle between α -C β -C β The angle between

[0039] 3.4) Calculate C based on the atomic three-dimensional coordinate information β -C β 、Tip-Tip、C α -Tip, Tip-C α The Euclidean distance of atom pairs is used to construct the spatial relationship matrix between residues;

[0040] 3.5) Use Foldseek to extract the 3D sequence features of amino acids and encode the spatial information in the protein structure with 20 letters: "A", "C", "D", "E", "F", "G", "H", "I", "K", "L", "M", "N", "P", "Q", "R", "S", "T", "V", "W", "Y";

[0041] 3.6) Calculate the The distance between them is used as the residue distance matrix feature.

[0042] The process of step 5) is as follows:

[0043] 5.1) Use a 3D convolutional network to process the voxelized features of the protein complex model. The process is as follows: the 3D voxelized features are processed through a 20×20 voxelization layer, and then processed through a 3D convolution layer with a convolution kernel size of 1×1, resulting in a sequence embedding feature of L×20×20×20×167, where L represents the length of the protein complex residues. By adjusting the step size of the convolution kernel, the 3D spatial features at the residue level are reduced to L×4×4×4×10. The GELU activation function is added after each convolution layer. After a series of pooling and horizontal and vertical striping processes, a 2D embedding vector of L×L×128 is finally obtained.

[0044] 5.2) A graph self-attention neural network (GAT) was used to process various features of the protein complex model. The process is as follows: the BLOSUM62 matrix, one-hot encoding features, main chain torsion angles, main chain bond lengths and angles, Rosetta residue intra-energy terms, secondary structure information, amino acid physical energy terms, 3di sequence features, protein structure sequence spectrum features, normalized distances between amino acids and sequence ends, and inter-residue orientation features were concatenated and used as GAT node features with a dimension of L×164. The residue distance matrix features were used as GAT edge features. After passing through a four-layer GAT network, each GAT layer passed through a normalization layer and a GELU activation function, and was pooled through a SAGPooling layer with a pooling rate of 0.5. Finally, after passing through a fully connected layer, a one-dimensional embedding vector of L×64 was obtained. The update formula for the current node feature of GAT is as follows:

[0045]

[0046] in, represents the new feature vector of node i after the GAT layer transformation, σ represents the activation function, represents the scaling factor, represents the attention coefficient, W k represents the weight matrix, represents the original feature vector of neighbor node j;

[0047] The calculation formula of the attention weight between two nodes in GAT is as follows:

[0048]

[0049] in, represents the attention weight of node i to neighbor node j under the kth attention head; represents the learnable weight vector unique to the kth attention head; W k represents the weight matrix associated with the kth attention head; It means that the features of node i and node j are transformed by the weight matrix of the kth attention head and spliced ​​together as the input of the attention mechanism;

[0050] 5.3) Using an axial self-attention neural network to process multiple features of the protein complex model, the process is as follows: the inter-residue distance map, the contact area and direction between residues, the Euler angles, and the two-dimensional energy distribution matrix are concatenated to obtain a feature matrix of dimension 32×L×L. This matrix is ​​then transformed to 128×L×L through two-dimensional convolution. Each convolutional layer is followed by a GELU activation function and a normalization layer. The matrix is ​​then concatenated with the two-dimensional embedding vector obtained in step 5.1). After two-dimensional convolution, normalization, and dimensionality reduction, a two-dimensional vector of dimension 1×128×L×L is obtained.

[0051] 5.4) Expand the shape obtained by calculating the mask in step 5.2) to 1×L×1×1, multiply it point-by-point with the vector obtained in step 5.3), and input it into a two-dimensional convolutional network with a convolution kernel size of 3×3 for decoding, and finally obtain the global quality score of the protein complex model.

[0052] The beneficial effects of the present invention are mainly manifested in: using a deep fusion network to evaluate the quality of protein models, fully integrating multidimensional features such as evolution, physics and geometry, and accurately capturing the complex topological relationships of protein complex models. This method uses three-dimensional convolution, graph self-attention and axial self-attention mechanisms to achieve efficient fusion of different feature spaces, enhance the structural characterization capabilities of key areas, and significantly improve the accuracy and robustness of quality assessment. Based on the deep fusion strategy, this method can more accurately predict the global quality score of protein complex models, provide reliable theoretical support for downstream structure optimization, molecular docking and protein function research, and contribute to the in-depth development of bioinformatics and structural biology. BRIEF DESCRIPTION OF THE DRAWINGS

[0053] Figure 1 This is an overall flow chart of a global quality assessment method for protein complex models based on deep fusion networks. DETAILED DESCRIPTION

[0054] The present invention will be further described below with reference to the accompanying drawings.

[0055] Reference Figure 1 , a global quality assessment method for protein complex models based on deep fusion network, comprising the following steps:

[0056] 1) Extract evolutionary features based on protein models. The process is as follows:

[0057] 1.1) Protein sequence alignment using the BLOSUM62 matrix. The process is as follows:

[0058] 1.1.1) Perform multiple sequence alignment on protein sequences and count the substitution frequencies of each pair of amino acids;

[0059] 1.1.2) Convert the frequency of each pair of amino acid substitutions into probabilities, assuming that the probability of substitution is P(X,Y), where X and Y are two amino acid types;

[0060] 1.1.3) Apply the logarithmic function to transform the substitution probability of each pair of amino acids and calculate their logarithmic ratio:

[0061]

[0062] 1.1.4) Fill the calculated logarithmic ratio values ​​of all amino acid pairs into the matrix to obtain the final BLOSUM62 matrix;

[0063] 1.2) Protein structure sequence spectrum characteristics, the process is:

[0064] 1.2.1) Using each protein complex model as a benchmark, use the Foldseek template search tool to search for template structures with a similarity score > 0.5 in the PAcluster80 template library, obtaining a total of n sets of matching templates;

[0065] 1.2.2) Perform structural alignment of the searched template structures using US-align, inserting gap marks ("-") for residue positions that do not form valid atomic correspondences;

[0066] 1.2.3) Based on C β Atomic spatial coordinate information is used to calculate the Euclidean distances of all valid residue pairs between the target structure and each template structure. For residue sites with gaps ("-"), the distance calculation of the position is automatically skipped, and finally a distance distribution matrix containing n+1 groups of data is generated;

[0067] 1.2.4) The distance range is divided into 36 equally spaced intervals (step length ), statistically analyze the probability distribution of the same position in each distance spectrum falling into each interval, and select the distance value corresponding to the probability peak as the structural sequence spectrum characteristic parameter of the site;

[0068] 1.3) Calculate the normalized distance between the amino acid and the end of the sequence;

[0069] 2) Extract physical features based on the protein model. The process is as follows:

[0070] 2.1) One-hot encoding is used to vectorize the amino acid sequence of the protein complex, generating a matrix containing 20 independent binary feature bits, each corresponding to a standard amino acid type;

[0071] 2.2) Extract the main chain torsion angle ( and ψ), and calculate the main chain C α -C, NCα , C-N bond length and C α -

[0072] CN、NC α -C, CNC α bond angle;

[0073] 2.3) Use pyRosetta to extract intra-residue energy features, including: intra-residue conformational preferences

[0074] (p_aa_pp), stereochemical constraints (rama_prepreo, omega), side chain conformational free energy

[0075] (fa_dun);

[0076] 2.4) Use pyRosetta to extract the interaction energy between residues, including van der Waals attraction (fa_atr) and repulsion (fa_rep), generalized solvation effect (fa_sol), and Lennard-Jones spherical core weighted potential.

[0077] (lk_ball_wtd), charge interaction (fa_elec), main chain-side chain hydrogen bond network (hbond_bb_sc)

[0078] and inter-side chain hydrogen bond synergy (hbond_sc);

[0079] 2.5) Analyze the secondary structure characteristics of protein complexes based on the DSSP algorithm;

[0080] 2.6) Construct a two-dimensional energy distribution matrix based on the distance between residues (C β -C β ) and relative orientation angle as coordinate axes, and the normalized energy distribution heat map is generated by Gaussian kernel density estimation;

[0081] 2.7) Extract the physicochemical properties of amino acids, including isoelectric point, polarity index, acidity / alkalinity, hydrogen bond donor / acceptor capacity, hydrophilicity / hydrophobicity, topological polar surface area, and secondary structure tendency;

[0082] 2.8) Use the Voronota tool to extract the contact area and orientation features between residues of the protein model based on the Voronoi tessellation;

[0083] 3) Extract geometric features based on the protein model. The process is as follows:

[0084] 3.1) Calculate the spatial voxel features of the protein structure: First, define a local coordinate system based on the positions of the main chain atoms with the current residue as the center, where the axis direction is determined by the vector pointing to the N atom of the current residue, the starting point is the C atom of the previous residue, and for the first residue, the axis direction points to the C atom of the current residue. α Atom. Calculate three orthogonal basis vectors by vector product and ensure that the relative position in the local coordinate system is fixed so that it is always in the axis direction and in the axis plane. Then, select the fixed coordinate point As the origin, construct a coordinate system with a side length of cube, and further divided into 20×20×20 The number of atoms in each voxel is counted and its center coordinates are calculated. At the same time, a value is assigned to each voxel to characterize the spatial density and distribution characteristics of the area.

[0085] 3.2) Calculate the Euler angle, using C α -C β The vector is the rotation axis, and its value relative to the reference coordinate system (NC α -C β α, β, and γ rotation angles in the direction of the plane normal;

[0086] 3.3) Calculate the positional relationship between residues, including: C α -C β -C β -C α The dihedral angle between NC α -C β -C β The dihedral angle between α -C β -C β The angle between

[0087] 3.4) Calculate C based on the atomic three-dimensional coordinate information β -C β 、Tip-Tip、C α -Tip, Tip-C α The Euclidean distance of atom pairs is used to construct the spatial relationship matrix between residues;

[0088] 3.5) Use Foldseek to extract the 3D sequence features of amino acids and encode the spatial information in the protein structure with 20 letters: "A", "C", "D", "E", "F", "G", "H",

[0089] "I", "K", "L", "M", "N", "P", "Q", "R", "S", "T", "V", "W", "Y";

[0090] 3.6) Calculate the The distance between them is used as the residue distance matrix feature;

[0091] 4) Producing label data: The label value is the template modeling score TM-score obtained by comparing each structure with the corresponding natural structure. The label value ranges from 0 to 1, where the label value corresponding to the natural structure is 1;

[0092] 5) After a series of preprocessing such as concatenation and dimensionality processing, the extracted features are input into the deep fusion network. The Adam optimizer (learning rate 0.001) is used to optimize the network weight parameters and the mean square error (MSE) loss function for training. Finally, the global score of the protein structure is calculated by the decoding module of the deep fusion network to further evaluate the global quality of the protein model. The process is as follows:

[0093] 5.1) Use a 3D convolutional network to process the voxelized features of the protein complex model. The process is as follows: the 3D voxelized features are processed through a 20×20 voxelization layer, followed by a 3D convolutional layer with a convolution kernel size of 1×1, resulting in a sequence embedding feature of L×20×20×20×167. By adjusting the stride of the convolution kernel, the 3D features at the residue level are reduced to L×4×4×4×10. The GELU activation function is added after each convolution layer. After a series of pooling and horizontal and vertical striping, a 2D embedding vector of L×L×128 is finally obtained.

[0094] 5.2) Use the GAT network to process various features of the protein complex model. The process is as follows: BLOSUM62 matrix, one-hot encoding features, main chain torsion angle, main chain bond length and angle, Rosetta residue intra-energy term, secondary structure information, physical energy term of amino acids, 3di sequence features, protein structure sequence spectrum features, normalized distance between amino acids and sequence ends, and inter-residue orientation features are spliced ​​and used as GAT node features (dimension is L×164), and the residue distance matrix features are used as GAT edge features. After passing through the four-layer GAT network, each layer of GAT passes through a normalization layer and GELU activation function, and is pooled through a SAGPooling layer with a pooling rate of 0.5. Finally, after passing through the fully connected layer, a one-dimensional embedding vector of L×64 is obtained. Among them, the current node of GAT is

[0095] The update formula of the feature is as follows:

[0096]

[0097] in, represents the new feature vector of node i after the GAT layer transformation, σ represents the activation function,

[0098] represents the scaling factor, represents the attention coefficient, W k represents the weight matrix, represents the original feature vector of neighbor node j;

[0099] The calculation formula of the attention weight between two nodes in GAT is as follows:

[0100]

[0101] in, represents the attention weight of node i to neighbor node j under the kth attention head;

[0102] represents the learnable weight vector unique to the kth attention head; W k represents the weight matrix associated with the kth attention head; It means that the features of node i and node j are transformed by the weight matrix of the kth attention head and spliced ​​together as the input of the attention mechanism;

[0103] 5.3) Use the axial self-attention neural network to process multiple features of the protein complex model. The process is as follows:

[0104] The distance map between residues, the contact area and direction between residues, the Euler angle, and the two-dimensional energy distribution matrix are concatenated to obtain a feature matrix with a dimension of 32×L×L. The dimension is transformed to 128×L×L through two-dimensional convolution. After each convolution layer, the GELU activation function and the normalization layer are passed, and then concatenated with the two-dimensional embedding vector obtained in step 6.1). After two-dimensional convolution, normalization, and dimensionality reduction, the dimension is obtained.

[0105] A two-dimensional vector of size 1×128×L×L;

[0106] 5.4) Expand the shape obtained by calculating the mask in step 6.2) to 1×L×1×1, multiply it point-by-point with the vector obtained in step 6.3), and input it into a two-dimensional convolutional network with a convolution kernel size of 3×3 for decoding, ultimately obtaining the global quality score of the protein complex model.

[0107] This example uses the natural protein structure 7QIH with an amino acid sequence length of 398 as an example. A global quality assessment method for protein complex models based on a deep fusion network includes the following steps:

[0108] 1) Extract evolutionary features based on protein models. The process is as follows:

[0109] 1.1) Protein sequence alignment substitution matrix BLOSUM62 matrix;

[0110] 1.2) Protein structure sequence spectrum characteristics, the process is;

[0111] 1.3) Calculate the normalized distance between the amino acid and the end of the sequence;

[0112] 2) Extract physical features based on the protein model. The process is as follows:

[0113] 2.1) One-hot encoding is used to vectorize the amino acid sequence of the protein complex;

[0114] 2.2) Extract the main chain torsion angle ( and ψ), and calculate the main chain C α -C, NC α , C-N bond length and C α -

[0115] CN、NC α -C, CNC α bond angle;

[0116] 2.3) Use pyRosetta to extract intra-residue energy features, including: intra-residue conformational preferences

[0117] (p_aa_pp), stereochemical constraints (rama_prepreo, omega), side chain conformational free energy

[0118] (fa_dun);

[0119] 2.4) Use pyRosetta to extract the interaction energy between residues, including van der Waals attraction (fa_atr) and repulsion (fa_rep), generalized solvation effect (fa_sol), and Lennard-Jones spherical core weighted potential.

[0120] (lk_ball_wtd), charge interaction (fa_elec), main chain-side chain hydrogen bond network (hbond_bb_sc)

[0121] and inter-side chain hydrogen bond synergy (hbond_sc);

[0122] 2.5) Analyze the secondary structure characteristics of protein complexes based on the DSSP algorithm;

[0123] 2.6) Construct a two-dimensional energy distribution matrix based on the distance between residues (C β -C β ) and relative orientation angle as coordinate axes, and the normalized energy distribution heat map is generated by Gaussian kernel density estimation;

[0124] 2.7) Extract the physicochemical properties of amino acids, including isoelectric point, polarity index, acidity / alkalinity, hydrogen bond donor / acceptor capacity, hydrophilicity / hydrophobicity, topological polar surface area, and secondary structure tendency;

[0125] 2.8) Use the Voronota tool to extract the contact area and orientation features between residues of the protein model based on the Voronoi tessellation;

[0126] 3) Extract geometric features based on the protein model. The process is as follows:

[0127] 3.1) Calculate the spatial voxel features of protein structure;

[0128] 3.2) Calculate the Euler angle, using C α -C β The vector is the rotation axis, and its value relative to the reference coordinate system (NC α -C β α, β, and γ rotation angles in the direction of the plane normal;

[0129] 3.3) Calculate the positional relationship between residues, including: C α -C β -C β -C α The dihedral angle between NC α -C β -C β The dihedral angle between α -C β -C β The angle between

[0130] 3.4) Calculate C based on the atomic three-dimensional coordinate information β -C β 、Tip-Tip、C α -Tip, Tip-C α The Euclidean distance of atom pairs is used to construct the spatial relationship matrix between residues;

[0131] 3.5) Use Foldseek to extract the 3D sequence features of amino acids and encode the spatial information in the protein structure with 20 letters: "A", "C", "D", "E", "F", "G", "H",

[0132] "I", "K", "L", "M", "N", "P", "Q", "R", "S", "T", "V", "W", "Y";

[0133] 3.6) Calculate the The distance between them is used as the residue distance matrix feature;

[0134] 4) Producing label data: The label value is the template modeling score TM-score obtained by comparing each structure with the corresponding natural structure. The label value ranges from 0 to 1, where the label value corresponding to the natural structure is 1;

[0135] 5) After a series of preprocessing such as concatenation and dimensionality processing, the extracted features are input into the deep fusion network. The Adam optimizer (learning rate 0.001) is used to optimize the network weight parameters and the mean square error (MSE) loss function for training. Finally, the global score of the protein structure is calculated by the decoding module of the deep fusion network to further evaluate the global quality of the protein model. The process is as follows:

[0136] 5.1) Using a 3D convolutional network to process the voxelized features of the protein complex model, the 3D voxelized features are processed through a 20×20 voxelization layer, followed by a 3D convolutional layer with a 1×1 kernel size, resulting in a sequence embedding feature of 398×20×20×20×167. By adjusting the stride of the convolution kernel, the 3D features at the residue level are reduced to 398×4×4×4×10. After each convolution layer, GELU activation is added. After a series of pooling and horizontal and vertical striping operations, a 2D embedding vector of 398×398×128 is finally obtained.

[0137] 5.2) A graph-based attention neural network (GAT) was used to process various features of the protein complex model. The process was as follows: the BLOSUM62 matrix, one-hot encoding features, main-chain torsion angles, main-chain bond lengths and angles, Rosetta residue intra-energy terms, secondary structure information, amino acid physical energy terms, 3-di sequence features, protein structure sequence spectrum features, normalized distances between amino acids and sequence ends, and inter-residue orientation features were concatenated and used as GAT node features (dimension 398 × 164). The residue distance matrix features were used as GAT edge features. The model was passed through a four-layer GAT network. Each GAT layer was followed by a normalization layer and a GELU activation function, and then pooled through a SAGPooling layer with a pooling rate of 0.5. Finally, a fully connected layer was used to obtain a 398 × 64 one-dimensional embedding vector.

[0138] 5.3) Use the axial self-attention neural network to process multiple features of the protein complex model. The process is as follows:

[0139] The inter-residue distance map, the contact area and direction between residues, the Euler angles, and the two-dimensional energy distribution matrix are concatenated to obtain a feature matrix of dimension 32 × 398 × 398. This matrix is ​​then transformed to 128 × 398 × 398 through two-dimensional convolution. Each convolutional layer is followed by a GELU activation function and a normalization layer. This matrix is ​​then concatenated with the two-dimensional embedding vector obtained in step 5.1. After two-dimensional convolution, normalization, and dimensionality reduction, a two-dimensional vector of dimension 1 × 128 × 398 × 398 is obtained.

[0140] 5.4) Expand the shape obtained by calculating the mask in step 5.2) to 1×398×1×1, multiply it point by point with the vector obtained in step 5.3), and input it into a two-dimensional convolutional network with a convolution kernel size of 3×3 for decoding, and finally obtain the global quality score of the protein complex model 7QIH.

[0141] Taking the natural structure 7QIH of a protein with an amino acid sequence length of 398 as an example, its true global quality score is 1. Using the above method, the global quality score of the protein complex model obtained is 0.9798, and the error of the result is very small.

[0142] The above description is about an embodiment of the present invention showing good effects. Obviously, the present invention is not only suitable for the above embodiment, but also can be implemented with various changes without departing from the basic spirit of the present invention and without exceeding the content involved in the essential content of the present invention.

Claims

1. A global quality assessment method for protein complex models based on deep fusion networks, characterized by: First, specific thresholds were set for maximum sequence redundancy, resolution, and residue length in the protein database (RCSB PDB) to screen out natural protein structures. For each protein structure, the HDock and xTrimoMultimer docking algorithms were used in combination with a chain orientation perturbation strategy to generate a perturbed structural model. Then, evolutionary, physical, and geometric features were extracted based on each protein structure model, and a deep fusion network was constructed that integrated a deep convolutional neural network, an axial attention mechanism network, and a graph attention network to fully capture the local and global topological information of the protein complex. Finally, the global quality score was predicted by the decoding module.

2. The global quality assessment method of protein complex model based on deep fusion network according to claim 1, characterized in that The method comprises the following steps: 1) Extract evolutionary features based on protein models. The process is as follows: 1.1) Protein sequence alignment substitution matrix BLOSUM62 matrix; 1.2) Protein structure sequence spectrum characteristics; 1.3) Calculate the normalized distance between the amino acid and the end of the sequence; 2) Extract physical features based on protein models; 3) Extracting geometric features based on protein models; 4) Producing label data: The label value is the template modeling score TM-score obtained by comparing each structure with the corresponding natural structure. The label value ranges from 0 to 1, where the label value corresponding to the natural structure is 1; 5) After concatenating and dimensionally processing the extracted features, they are input into the deep fusion network. The Adam optimizer is used to optimize the network weight parameters and the mean square error (MSE) loss function for training. Finally, the global score of the protein structure is calculated by the decoding module of the deep fusion network to further evaluate the global quality of the protein model.

3. The global quality assessment method of protein complex model based on deep fusion network according to claim 2, characterized in that: The process of 1.1) is as follows: 1.1.1) Perform multiple sequence alignment on protein sequences and count the substitution frequencies of each pair of amino acids; 1.1.2) Convert the frequency of each pair of amino acid substitutions into probabilities, assuming that the probability of substitution is P(X,Y), where X and Y are two amino acid types; 1.1.3) Apply the logarithmic function to transform the substitution probability of each pair of amino acids and calculate their logarithmic ratio: 1.1.4) Fill the calculated logarithmic ratio values ​​of all amino acid pairs into the matrix to obtain the final BLOSUM62 matrix.

4. The global quality assessment method of protein complex model based on deep fusion network according to claim 2, characterized in that The process of 1.2) is as follows: 1.2.1) Using each protein complex model as a benchmark, use the Foldseek template search tool to search for template structures with a similarity score > 0.5 in the PAcluster80 template library, obtaining a total of n sets of matching templates; 1.2.2) Perform structural alignment of the searched template structures using US-align, inserting gap marks ("-") for residue positions that do not form valid atomic correspondences; 1.2.3) Based on C β Atomic spatial coordinate information is used to calculate the Euclidean distances of all valid residue pairs between the target structure and each template structure. For residue sites with gaps ("-"), the distance calculation of the position is automatically skipped, and finally a distance distribution matrix containing n+1 groups of data is generated; 1.2.4) The distance range is divided into 36 equally spaced intervals with a step size of The probability distribution of the same position in each distance map falling into each interval is counted, and the distance value corresponding to the probability peak is selected as the structural sequence spectrum characteristic parameter of the site.

5. The global quality assessment method of protein complex model based on deep fusion network according to claim 2, characterized in that: The process of 2) is as follows: 2.1) One-hot encoding is used to vectorize the amino acid sequence of the protein complex, generating a matrix containing 20 independent binary feature bits, each corresponding to a standard amino acid type; 2.2) Extract the main chain torsion angle, and ψ, and calculate the main chain C α -C, NC α , C-N bond length and C α -CN、NC α -C, CNC α bond angle; 2.3) Use pyRosetta to extract residue intra-energy features, including intra-residue conformational preferences (p_aa_pp), stereochemical constraints (rama_prepreo, omega), and side chain conformational free energy (fa_dun); 2.4) Use pyRosetta to extract the interaction energies between residues, including van der Waals attraction (fa_atr) and repulsion (fa_rep), generalized solvation effect (fa_sol), Lennard-Jones spherical core weighted potential (lk_ball_wtd), charge interaction (fa_elec), main chain-side chain hydrogen bond network (hbond_bb_sc) and side chain hydrogen bond cooperativity (hbond_sc); 2.5) Analyze the secondary structure characteristics of protein complexes based on the DSSP algorithm; 2.6) Construct a two-dimensional energy distribution matrix based on the distance between residues (C β -C β ) and relative orientation angle as coordinate axes, and the normalized energy distribution heat map is generated by Gaussian kernel density estimation; 2.7) Extract the physicochemical properties of amino acids, including isoelectric point, polarity index, acidity / alkalinity, hydrogen bond donor / acceptor capacity, hydrophilicity / hydrophobicity, topological polar surface area, and secondary structure tendency; 2.8) Use the Voronota tool to extract the contact area and direction features between residues of the protein model based on Voronoi tessellation.

6. The global quality assessment method of protein complex model based on deep fusion network according to claim 2, characterized in that: The process of 3) is as follows: 3.1) Calculate the spatial voxel features of the protein structure: First, define a local coordinate system based on the positions of the main chain atoms with the current residue as the center, where the axis direction is determined by the vector pointing to the N atom of the current residue, the starting point is the C atom of the previous residue, and for the first residue, the axis direction points to the C atom of the current residue. α Atoms; calculate three orthogonal basis vectors by vector product and ensure that the relative position in the local coordinate system is fixed so that it is always in the axis direction and in the axis plane; then, select a fixed coordinate point As the origin, construct a coordinate system with a side length of cube, and further divided into 20×20×20 The number of atoms in each voxel is counted and its center coordinates are calculated. At the same time, a value is assigned to each voxel to characterize the spatial density and distribution characteristics of the area. 3.2) Calculate the Euler angle, using C α -C β The vector is the rotation axis, and its value is calculated relative to the reference coordinate system, i.e., NC α -C β α, β, and γ rotation angles in the direction of the plane normal; 3.3) Calculate the positional relationship between residues, including: C α -C β -C β -C α The dihedral angle between NC α -C β -C β The dihedral angle between α -C β -C β The angle between 3.4) Calculate C based on the atomic three-dimensional coordinate information β -C β 、Tip-Tip、C α -Tip, Tip-C α The Euclidean distance of atom pairs is used to construct the spatial relationship matrix between residues; 3.5) Use Foldseek to extract the 3D sequence features of amino acids and encode the spatial information in the protein structure with 20 letters: "A", "C", "D", "E", "F", "G", "H", "I", "K", "L", "M", "N", "P", "Q", "R", "S", "T", "V", "W", "Y"; 3.6) Calculate the The distance between them is used as the residue distance matrix feature.

7. The method for global quality assessment of protein complex models based on deep fusion network according to claim 2, characterized in that: The process of step 5) is as follows: 5.1) Use a 3D convolutional network to process the voxelized features of the protein complex model. The process is as follows: the 3D voxelized features are processed through a 20×20 voxelization layer, and then processed through a 3D convolution layer with a convolution kernel size of 1×1, resulting in a sequence embedding feature of L×20×20×20×167, where L represents the length of the protein complex residues. By adjusting the step size of the convolution kernel, the 3D spatial features at the residue level are reduced to L×4×4×4×10. The GELU activation function is added after each convolution layer. After a series of pooling and horizontal and vertical striping processes, a 2D embedding vector of L×L×128 is finally obtained. 5.2) Use the graph self-attention neural network (GAT) to process various features of the protein complex model. The process is as follows: the BLOSUM62 matrix, one-hot encoding features, main chain torsion angles, main chain bond lengths and angles, Rosetta residue intra-energy terms, secondary structure information, amino acid physical energy terms, 3di sequence features, protein structure sequence spectrum features, normalized distances between amino acids and sequence ends, and inter-residue orientation features are spliced ​​and used as GAT node features with a dimension of L×164. The residue distance matrix features are used as GAT edge features. After passing through a four-layer GAT network, each GAT layer passes through a normalization layer and a GELU activation function, and is pooled through a SAGPooling layer with a pooling rate of 0.

5. Finally, after passing through a fully connected layer, a one-dimensional embedding vector of L×64 is obtained, where The update formula of the current node feature of GAT is as follows: in, represents the new feature vector of node i after the GAT layer transformation, σ represents the activation function, represents the scaling factor, represents the attention coefficient, W k represents the weight matrix, represents the original feature vector of neighbor node j; The calculation formula of the attention weight between two nodes in GAT is as follows: in, represents the attention weight of node i to neighbor node j under the kth attention head; represents the learnable weight vector unique to the kth attention head; W k represents the weight matrix associated with the kth attention head; It means that the features of node i and node j are transformed by the weight matrix of the kth attention head and spliced ​​together as the input of the attention mechanism; 5.3) Using an axial self-attention neural network to process multiple features of the protein complex model, the process is as follows: the inter-residue distance map, the contact area and direction between residues, the Euler angles, and the two-dimensional energy distribution matrix are concatenated to obtain a feature matrix of dimension 32×L×L. This matrix is ​​then transformed to 128×L×L through two-dimensional convolution. Each convolutional layer is followed by a GELU activation function and a normalization layer. The matrix is ​​then concatenated with the two-dimensional embedding vector obtained in step 5.1). After two-dimensional convolution, normalization, and dimensionality reduction, a two-dimensional vector of dimension 1×128×L×L is obtained. 5.4) Expand the shape obtained by calculating the mask in step 5.2) to 1×L×1×1, multiply it point-by-point with the vector obtained in step 5.3), and input it into a two-dimensional convolutional network with a convolution kernel size of 3×3 for decoding, and finally obtain the global quality score of the protein complex model.

Images