Missense mutation function effect prediction model construction method and prediction method

Abstract

The invention discloses a missense mutation function effect prediction model construction method and a prediction method, and relates to the technical field of biological information. According to the method, on the basis of pathogenic and benign missense mutation data, various molecular phenotypic characteristics are enriched and annotated through statistical indexes, and a machine learning model for mutation function effect prediction is constructed. The model fuses feature information of a plurality of molecular levels, comprises comprehensive evaluation indexes of a residual base level, a mutation level and a protein level, and is trained and optimized based on a random forest algorithm. The method not only can accurately predict pathogenic mutation, but also shows excellent performance in recognition of cancer-driven missense mutation, has good generalization ability and biological interpretability, and can provide solid technical support for pathogenicity and carcinogenicity recognition of missense mutation and molecular mechanism research of missense mutation.

Description

Technical Field

[0001] The present invention relates to a method for constructing a missense mutation functional effect prediction model and a prediction method, belonging to the field of bioinformatics technology. Background Art

[0002] In recent years, the development of high-throughput sequencing technology has revealed a large number of missense mutations (hereinafter referred to as mutations) in the human genome, but only some of them play a key role in the occurrence and development of diseases or cancers. Therefore, accurately identifying pathogenic mutations and cancer driver mutations with functional impacts has become an important scientific problem to be solved urgently. Although functional verification experiments can provide direct and reliable evidence, due to high costs and long cycles, it is difficult to be applied on a large scale to screen all potential mutations.

[0003] For this reason, a variety of computational prediction methods and tools have been developed in recent years to evaluate the functions and phenotypic impacts of missense mutations. However, most of the existing methods are only constructed for pathogenic mutations and it is difficult to simultaneously take into account the identification of cancer driver mutations, especially showing limitations in driver judgment. In addition, many methods rely too much on fitness-based scoring systems and lack biological interpretability, making it difficult to reveal potential molecular mechanisms.

[0004] Conceptually, the "functional effect" of a mutation usually refers to whether this mutation has a significant impact on protein structure, function or interaction. A pathogenic mutation means that this molecular effect ultimately leads to abnormal body phenotypes or disease manifestations; while a cancer driver mutation generates a selective advantage through obtaining functional changes during the evolution of cancer cells. Therefore, the functional effect is the common basis for pathogenicity and drivability, and there is significant overlap between the two types of mutations at the molecular level. Based on this commonality, constructing a unified mutation effect prediction model has become the key direction to improve prediction accuracy and generalization ability. Summary of the Invention

[0005] The purpose of the present invention is to overcome the deficiencies in the prior art and provide a method for constructing a missense mutation functional effect prediction model and a prediction method, which can not only accurately predict pathogenic mutations, but also show excellent performance in the identification of cancer driver missense mutations, have good generalization ability and biological interpretability, and can provide solid technical support for the identification of the pathogenicity and carcinogenicity of missense mutations and the research on their molecular mechanisms.

[0006] To achieve the above object, the present invention is implemented by adopting the following technical solutions:

[0007] On the one hand, the present invention provides a method for constructing a missense mutation functional effect prediction model, including:

[0008] Integrating missense mutation data from different sources to construct a training data set, where the missense mutations include pathogenic missense mutations and benign missense mutations;

[0009] For each missense mutation data in the training dataset, perform enrichment analysis on molecular phenotype features and protein-level features respectively to obtain the molecular phenotype feature evaluation value and the comprehensive evaluation value of the major protein-level features, and calculate the predicted values of the protein stability change and binding affinity change caused by the missense mutation;

[0010] According to the molecular phenotype feature evaluation value and the comprehensive evaluation value of the major protein-level features, perform molecular phenotype feature annotation and major protein feature annotation on each missense mutation in the training dataset to obtain the molecular phenotype feature annotation result and the major protein-level feature annotation result respectively;

[0011] Use the molecular phenotype feature annotation result, the major protein-level feature annotation result of each missense mutation in the training dataset, and the predicted values of the protein stability change and binding affinity change caused by the missense mutation to form a multi-dimensional feature vector;

[0012] Construct a random forest classification model;

[0013] Use the multi-dimensional feature vector as the input to train the random forest classification model to obtain a missense mutation functional effect prediction model.

[0014] Furthermore, the missense mutation data includes an identifier, a missense mutation position, a reference residue, and a mutant residue.

[0015] Furthermore, after constructing the training dataset, it also includes missense mutation standardization and protein structure mapping for the training dataset;

[0016] The missense mutation standardization is used to unify the missense mutation data with the standard protein name of the unified protein resource database as the identifier, specifically including:

[0017] When the identifier of the missense mutation is a gene symbol, convert the gene symbol where the missense mutation is located into the standard protein name of the unified protein resource database, and check whether the missense mutation position and the reference residue are consistent with the residues at the same position in the standard protein sequence of the unified protein resource database. If they are consistent, retain them; if not, delete the missense mutation;

[0018] When the identifier of the missense mutation is the standard protein name of the unified protein resource database, check whether the missense mutation position and the reference residue are consistent with the residues at the same position in the standard protein sequence of the unified protein resource database. If they are consistent, retain them; if not, delete the missense mutation;

[0019] When the identifier of a missense mutation is a transcript in a reference sequence database, a transcript in a genomic database, or a protein in a genomic database, convert the transcript in the reference sequence database, the transcript in the genomic database, or the protein in the genomic database where the missense mutation is located into the isoform protein name in the unified protein resource database, and check whether the position of the missense mutation and the reference residue are consistent with the residues at the same position in the isoform protein sequence of the unified protein resource database. If they are consistent, retain them; if not, delete the mutation.

[0020] Perform a pairwise sequence alignment between the standard protein sequence and the isoform protein sequence in the unified protein resource database, map the mutation position originally based on the isoform protein sequence in the unified protein resource database to the mutation position based on the standard protein sequence in the unified protein resource database, and check whether the reference residue of the missense mutation is consistent with the residue at the corresponding missense mutation position in the standard protein sequence of the unified protein resource database. If they are consistent, retain them; if not, delete the missense mutation.

[0021] The retained missense mutations are the standardized missense mutations.

[0022] The protein structure mapping includes:

[0023] Map the missense mutation to the AlphaFold2 structure sequence according to the position and reference residue of the standardized missense mutation in the standard protein sequence of the unified protein resource database. Among them, if the standard protein sequence of the unified protein resource database of the protein where the missense mutation is located is the same as the AlphaFold2 structure sequence, map it directly; if the standard protein sequence of the unified protein resource database of the protein where the missense mutation is located is different from the AlphaFold2 structure sequence, perform a pairwise sequence comparison between the standard protein sequence of the unified protein resource database and the AlphaFold2 structure sequence to complete the mapping.

[0024] Furthermore, the enrichment analysis of the molecular phenotype characteristics of each missense mutation data in the training dataset to obtain the molecular phenotype characteristic evaluation value includes:

[0025] Calculate the number of pathogenic missense mutations and benign missense mutations that conform to the molecular phenotype characteristics and the number of pathogenic missense mutations and benign missense mutations that do not conform to the molecular phenotype characteristics in the training dataset.

[0026] Take the number of pathogenic missense mutations and benign missense mutations that conform to the molecular phenotype characteristics and the number of pathogenic missense mutations and benign missense mutations that do not conform to the molecular phenotype characteristics as inputs, and calculate the odds ratio of pathogenic mutations and benign mutations in this molecular phenotype characteristic based on the two-sided Fisher's exact test method.

[0027] Take the base-10 logarithmic function of the odds ratio of the molecular phenotypic feature as the evaluation value of the molecular phenotypic feature.

[0028] Furthermore, the molecular phenotypic features are divided into residue-level features and mutation-level features;

[0029] The residue-level features include that the residue where the mutation occurs is a core residue, the residue where the mutation occurs is a surface residue, the secondary structure of the protein where the mutation occurs is a helix type, the secondary structure of the protein where the mutation occurs is a fold type, the secondary structure of the protein where the mutation occurs is a loop structure type, the mutation occurs at a post-translational modification site, the residue within a spatial distance of 4 Å from the post-translational modification site where the mutation occurs, the residue within a spatial distance of 8 Å from the post-translational modification site where the mutation occurs, the mutation occurs at a disordered region site, the mutation occurs at a short linear motif site, the mutation occurs at an allosteric site, the mutation occurs at a protein-protein binding site, the mutation occurs at a protein-small molecule ligand binding site, the mutation occurs at a protein-nucleic acid binding site, the mutation occurs at a protein domain site, the physicochemical property of the mutation reference residue belongs to aliphatic, the physicochemical property of the mutation reference residue belongs to aromatic, the physicochemical property of the mutation reference residue belongs to positively charged, the physicochemical property of the mutation reference residue belongs to negatively charged, the physicochemical property of the mutation reference residue belongs to neutral, the physicochemical property of the mutation reference residue belongs to special;

[0030] The mutation-level features include that the absolute value of the charge change of the amino acid residue before and after the mutation is 0, the absolute value of the charge change of the amino acid residue before and after the mutation is 1, the absolute value of the charge change of the amino acid residue before and after the mutation is 2, the absolute value of the volume change of the amino acid residue before and after the mutation is 0, the absolute value of the volume change of the amino acid residue before and after the mutation is 1, the absolute value of the volume change of the amino acid residue before and after the mutation is 2.

[0031] Furthermore, the physicochemical properties of the mutation reference residue are specifically classified based on 20 standard amino acids as follows:

[0032] Alanine, isoleucine, leucine, methionine, and valine belong to aliphatic; phenylalanine, tryptophan, and tyrosine belong to aromatic; histidine, lysine, and arginine belong to positively charged; aspartic acid and glutamic acid belong to negatively charged; asparagine, glutamine, serine, and threonine belong to neutral; cysteine, proline, and glycine belong to special;

[0033] The method for judging the absolute value of the charge change of the amino acid residue before and after the mutation includes:

[0034] Set charge parameters, specifically including: lysine and arginine are +1 positive charge, aspartic acid and glutamic acid are -1 negative charge, and the remaining amino acids are 0 neutral charge;

[0035] Calculate the absolute value of the charge change of the amino acid residue before and after mutation according to the charge parameter and the amino acid type of the amino acid residue before and after mutation.

[0036] The method for judging the absolute value of the volume change of the amino acid residue before and after mutation includes:

[0037] Set volume coding, specifically including: the volume coding of alanine, glycine, serine, cysteine, proline, threonine, aspartic acid and asparagine is 0, the volume coding of valine, histidine, glutamic acid and glutamine is 1, and the volume coding of isoleucine, leucine, methionine, lysine, arginine, phenylalanine, tryptophan and tyrosine is 2.

[0038] Calculate the absolute value of the volume change of the amino acid residue before and after mutation according to the volume coding and the amino acid type of the amino acid residue before and after mutation.

[0039] Further, the large-class protein level features include human protein function types, human protein subcellular localization information types, the standard protein range of the Kyoto Encyclopedia of Genes and Genomes pathways and the unified protein resource database corresponding to cancer signaling pathways, and the protein specific domain range.

[0040] Each large-class protein level feature includes multiple small-class protein level features, and the comprehensive evaluation value of the large-class protein level feature is obtained by summing the evaluation values of all small-class protein level features that the missense mutation in the protein conforms to in each large-class protein level feature.

[0041] Further, the enrichment analysis of the protein level features for each missense mutation in the training dataset to obtain the comprehensive evaluation value of the large-class protein level features includes:

[0042] Calculate the number of pathogenic missense mutations and benign missense mutations that conform to the small-class protein level features and the number of pathogenic missense mutations and benign missense mutations that do not conform to the small-class protein level features in the training dataset, and use them as inputs to calculate the odds ratio of each small-class protein level feature based on the two-sided Fisher's exact test method.

[0043] For any missense mutation, sum the logarithms of the odds ratios of all small-class protein level features that the protein where the missense mutation is located conforms to in each large-class protein level feature with base 10, and use it as the comprehensive evaluation value of the missense mutation in the corresponding large-class protein level feature.

[0044] Further, the calculation expression of the odds ratio of the missense mutation in the molecular phenotype feature is:

[0045] ;

[0046] Among them, represents the odds ratio of pathogenic missense mutations and benign missense mutations in a molecular phenotypic feature, represents the number of pathogenic missense mutations that conform to the molecular phenotypic feature, represents the number of pathogenic missense mutations that do not conform to the molecular phenotypic feature, represents the number of benign missense mutations that conform to the molecular phenotypic feature, represents the number of benign missense mutations that do not conform to the molecular phenotypic feature;

[0047] The calculation expression for the odds ratio of missense mutations in the small-class protein level feature is: ;

[0048] Among them, represents the odds ratio of pathogenic missense mutations and benign missense mutations in a small-class protein level feature, represents the number of pathogenic missense mutations that conform to the small-class protein level feature, represents the number of pathogenic missense mutations that do not conform to the small-class protein level feature, represents the number of benign missense mutations that conform to the small-class protein level feature, represents the number of benign missense mutations that do not conform to the small-class protein level feature;

[0049] The calculation expression for the comprehensive evaluation value of the large-class protein level feature is: ;

[0050] Among them, represents any one of the large-class protein level features, represents the comprehensive value of the odds ratio of pathogenic missense mutations and benign missense mutations in the large-class protein level feature, is the number of all small-class features that the protein where the mutation is located conforms to in the large-class protein feature, is the th odds ratio of the small-class feature.

[0051] Furthermore, based on the evaluation value of the molecular phenotypic feature and the comprehensive evaluation value of the large-class protein level feature, molecular phenotypic feature annotation and large-class protein feature annotation are performed on each missense mutation in the training data set to obtain a molecular phenotypic feature annotation result and a large-class protein level feature annotation result respectively, specifically including:

[0052] For any missense mutation in the training dataset, label the molecular phenotype features that the missense mutation conforms to as the evaluation values of the corresponding molecular phenotype features, and label the molecular phenotype features that do not conform to as 0, to obtain the molecular phenotype annotation result; label the small-class protein level features that the missense mutation conforms to as the evaluation values of the corresponding small-class protein level features, and label the small-class protein level features that do not conform to as 0, to obtain the annotation result of the small-class protein level features, and sum up the annotation results of the small-class protein level features corresponding to obtain the annotation result of the large-class protein level features.

[0053] On the other hand, the present invention also provides a method for predicting the functional effect of a missense mutation, including:

[0054] Obtain the position, reference residue, and mutant residue of the missense mutation to be tested in the standard sequence of the unified protein resource database;

[0055] According to the position, reference residue, and mutant residue of the missense mutation to be tested in the standard sequence of the unified protein resource database, determine whether the missense mutation to be tested conforms to the molecular phenotype features. If it does not conform to a certain molecular phenotype feature, the annotation result of this molecular phenotype feature is 0. If it conforms to a certain molecular phenotype feature, the annotation result of this molecular phenotype feature is the evaluation value of this molecular phenotype feature, to obtain the molecular phenotype feature judgment result;

[0056] According to the position, reference residue, and mutant residue of the missense mutation to be tested in the standard sequence of the unified protein resource database, for each large-class protein level feature, sum up the evaluation values of all small-class protein level features that the protein where the missense mutation to be tested is located conforms to in each large-class protein level feature to obtain the comprehensive evaluation value of this large-class protein level feature, and use the comprehensive evaluation values of all large-class protein level features as the protein level feature judgment result;

[0057] Calculate the predicted values of the protein stability change and binding affinity change caused by the missense mutation to be tested;

[0058] Use the molecular phenotype feature judgment result, protein level feature judgment result, and the predicted values of protein stability change and binding affinity change as inputs, and use the missense mutation functional effect prediction model obtained based on any one of the above-mentioned missense mutation functional effect prediction model construction methods to predict the functional effect of the missense mutation to be tested, and output to obtain the functional effect classification result of the missense mutation to be tested.

[0059] Further, the functional effect classification result includes: when the classification label is 1, it means that the missense mutation to be tested is a missense mutation with functional impact; when the classification label is 0, it means that the missense mutation to be tested is a missense mutation without functional impact;

[0060] Among them, the missense mutations with functional impacts include pathogenic missense mutations or cancer driver missense mutations; the missense mutations without functional impacts include benign missense mutations or passenger missense mutations.

[0061] Compared with the prior art, the beneficial effects achieved by the present invention are as follows:

[0062] The present invention combines the evaluation values of 27 molecular phenotype features, the comprehensive evaluation values of 5 types of protein level features, and the predicted values of the changes in protein stability and binding affinity caused by mutations to form 34 features, constructs a random forest classification model to obtain a missense mutation functional effect prediction model, and predicts the functional effects of missense mutations through the missense mutation functional effect prediction model. It can not only accurately predict pathogenic mutations, but also shows excellent performance in the identification of cancer driver missense mutations, has good generalization ability and biological interpretability, and can provide strong technical support for the identification of the pathogenicity and carcinogenicity of missense mutations and the research on their molecular mechanisms. It has broad application prospects in precision medicine fields such as drug target screening, clinical diagnosis, and personalized treatment. BRIEF DESCRIPTION OF THE DRAWINGS

[0063] Figure 1 It is a schematic flow chart of the construction method of a missense mutation functional effect prediction model based on molecular phenotype features in an embodiment of the present invention;

[0064] Figure 2 It is a schematic diagram of the analysis results of the log 10 OR Pathogenic / Benign of 27 molecular phenotype features in the construction method of a missense mutation functional effect prediction model based on molecular phenotype features in an embodiment of the present invention;

[0065] Figure 3 It is a schematic diagram comparing the test results of the ROC curves (Receiver Operating Characteristic curves) of the missense mutation functional effect prediction model in an embodiment of the present invention and other models on a pathogenic mutation independent test dataset;

[0066] Figure 4 It is a schematic diagram comparing the test results of the ROC curves of the missense mutation functional effect prediction model in an embodiment of the present invention and other models on a cancer driver mutation independent test dataset;

[0067] Figure 5 It is a schematic diagram of the changes in protein stability and binding affinity caused by pathogenic missense mutations (Pathogenic) and benign missense mutations (Benign) in an embodiment of the present invention. Among them, a is a schematic diagram of the change in protein stability, and b is a schematic diagram of the change in protein binding affinity. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0068] The present invention will be further described below in conjunction with the accompanying drawings. The following embodiments are only used to more clearly illustrate the technical solution of the present invention and should not be used to limit the protection scope of the present invention.

[0069] Embodiment 1:

[0070] As Figure 1 shown, the embodiment of the present invention provides a method for constructing a missense mutation functional effect prediction model based on molecular phenotypic characteristics, including the following steps:

[0071] First, obtain pathogenic / benign missense mutations from different sources:

[0072] In this embodiment, missense mutation data is collected from three authoritative data sources: the ClinVar database (Clinical Variation Database, hereinafter referred to as ClinVar) (May 2021), the UniProt Humsavar database (UniProt Humsavar Human Variation Database, hereinafter referred to as Humsavar) (February 2021), and gain-of-function and loss-of-function mutations sorted out by Itan Laboratory based on the December 2019 version of the HGMD database (Human Gene Mutation Database, hereinafter referred to as HGMD).

[0073] Preliminary screening is carried out on the collected ClinVar mutations and Humsavar mutations: only retain ClinVar missense mutations with clinical significance of "pathogenic / likely pathogenic" or "benign / likely benign", and at the same time exclude mutations with contradictory clinical significance; only retain Humsavar missense mutations labeled as "likely pathogenic or pathogenic" and "likely benign or benign".

[0074] Classify the "pathogenic / likely pathogenic" mutations in ClinVar, the "pathogenic / likely pathogenic" mutations in Humsavar, and all HGMD gain-of-function and loss-of-function mutations as pathogenic mutations; classify the "likely benign or benign" mutations in ClinVar and the "likely benign or benign" mutations in Humsavar as benign mutations.

[0075] Since missense mutation identifiers from different data sources are different. For example, the identifier of ClinVar mutation is the transcript of the reference sequence database (NCBI Reference Sequence Database, RefSeq), the identifier of Humsavar mutation is the standard protein of the Universal Protein Resource database (UniProt), and the identifiers of HGMD gain-of-function and loss-of-function mutations are the proteins of the genomic database (Ensembl). Therefore, it is necessary to standardize the mutation identifiers and their corresponding mutation positions. Given that the subsequent statistical annotations of this invention are all based on proteins, it is necessary to standardize the mutation data from different sources based on the UniProt standard protein. Since genes may produce multiple protein isoforms through mechanisms such as alternative splicing, RNA editing, or genomic rearrangement, among which the canonical sequence is the most representative and functionally important sequence among all naturally occurring isoforms. Therefore, all missense mutations are converted into mutations on the UniProt standard protein sequence (referred to as UniProt standard mutations). The specific method for standardizing missense mutations is as follows:

[0076] For mutations with the original identifier being the gene symbol, first use the identifier matching tool (ID Mapping) of the Universal Protein Resource database (UniProt) to convert the gene symbol into the UniProt standard protein name, and then check whether the provided mutation position and reference residue are consistent with the residue at the same position in the UniProt standard protein sequence. If they are consistent, retain the mutation; if not, delete the mutation.

[0077] For mutations with the original identifier being the UniProt standard protein name, directly check whether the provided mutation position and reference residue are consistent with the residue at the same position in the UniProt standard protein sequence. If they are consistent, retain the mutation; if not, delete the mutation.

[0078] For mutations with original identifiers of RefSeq transcripts, Ensembl transcripts, or Ensembl proteins, it is necessary to first convert the RefSeq transcript, Ensembl transcript, or Ensembl protein where the missense mutation is located into a UniProt isoform protein name using the Proteins API module provided by EBI (the protein-related data interface provided by EBI); secondly, check whether the provided mutation position and reference residue are consistent with the residue at the same position in the UniProt isoform sequence. If they are consistent, retain the mutation; if not, remove the mutation; finally, use Clustal Omega (a multiple sequence alignment tool) to perform a pairwise sequence alignment between the UniProt isoform protein sequence and the corresponding UniProt canonical protein sequence, convert the mutation position on the UniProt isoform sequence to the position on the UniProt canonical protein sequence, and check whether the reference residue of the mutation is consistent with the residue at the mutation position on the UniProt canonical protein sequence. If they are consistent, retain the mutation; if not, delete the mutation.

[0079] After completing the mutation standardization, observe whether there are identical UniProt canonical mutations that appear repeatedly among different data sources. If the same UniProt canonical mutation is classified into mutually exclusive categories in different sources (for example, this UniProt canonical mutation is a benign missense mutation in source 1 and a pathogenic missense mutation in source 2), it needs to be excluded; if the same UniProt canonical mutation is classified into the same category in different sources, retain this mutation in the dataset, and only retain the repeated UniProt canonical mutation once.

[0080] For the mutations in the dataset, according to subsequent statistical analyses, overly long proteins may introduce significant biases in the calculation of statistical metrics. Therefore, a statistical analysis of the sequence lengths of the proteins in the dataset was performed, and it was found that proteins with sequence lengths exceeding 3,000 only accounted for a small fraction of all proteins. Thus, proteins with sequence lengths exceeding 3,000 were removed, while retaining the vast majority of proteins and avoiding potential biases that may be introduced by overly long proteins in subsequent calculations.

[0081] Based on the above data processing, perform protein structure mapping on the missense mutations to obtain the missense mutations after protein structure mapping. Specifically:

[0082] In this embodiment, 18,531 three-dimensional structures of human proteins predicted by AlphaFold2 were downloaded from the AlphaFold database (the name of the protein structure prediction tool). To evaluate the confidence of the structure prediction, AlphaFold2 uses the pLDDT (predicted Local Distance Difference Test) score to describe the prediction accuracy of each amino acid residue in the structure. The value range of this score is from 0 to 100, and the higher the score, the more reliable the conformational prediction of the residue. Therefore, a screening criterion with the median of pLDDT of all residues in the predicted structure ≥ 70 was set, and finally 13,880 reliable protein structures were retained. Subsequently, according to the position of the mutation in the UniProt standard sequence and the reference residue, the mutation was mapped to the AlphaFold2 predicted structure. It was found that the AlphaFold2 structure sequences of most proteins were consistent with the UniProt standard sequence, and the mapping of protein sequences and structures could be directly achieved. However, there were still a small number of proteins with inconsistent sequences. For these proteins, Clustal Omega (a multiple sequence alignment tool) was used to perform a pairwise sequence alignment between the AlphaFold2 structure sequence and the UniProt standard sequence to achieve the structural mapping of the corresponding mutation.

[0083] So far, a total of 67,188 missense mutations successfully mapped to the structure have been obtained, including 29,844 pathogenic mutations and 37,344 benign mutations, which constitute the training dataset.

[0084] Then, based on the training dataset composed of pathogenic and benign missense mutations, enrichment analysis annotation of molecular phenotype characteristics was carried out using statistical indicators to obtain the evaluation value of molecular phenotype characteristics.

[0085] Molecular phenotype characteristics are divided into residue-level characteristics and mutation-level characteristics.

[0086] In this embodiment, a total of 21 residue-level characteristics were collected and calculated, as follows:

[0087] Residue exposure level: In the protein structure, according to the solvent accessible surface area (SASA), residues can be classified into two categories: core and surface. When the ratio of the SASA of a residue in the structure to its SASA in the extended tripeptide is less than 0.2, it is defined as a core residue; when it is greater than 0.2, it is a surface residue. The SASA of a residue in the structure can be calculated by DSSP (the name of the calculation tool), and its SASA in the extended tripeptide can be obtained from the literature.

[0088] Therefore, the first to second residue-level features are that the residue where the mutation occurs is a core residue and the residue where the mutation occurs is a surface residue, respectively.

[0089] Secondary structure: The secondary structure of a protein is identified by the DSSP program (name of the computational tool) and is classified into three categories: helix, strand, and loop.

[0090] Therefore, the third to fifth residue-level features are that the secondary structure of the protein where the mutation occurs is of the helix type, the secondary structure of the protein where the mutation occurs is of the strand type, and the secondary structure of the protein where the mutation occurs is of the loop type, respectively.

[0091] Post-Translational Modifications (PTM) sites and residues within a spatial distance of 4 Å and 8 Å: In this example, experimentally verified PTM sites were collected from two authoritative databases, PhosphoSitePlus and dbPTM, and all the collected PTM data were integrated, with duplicate PTM sites retained only once, resulting in 703,301 PTM site feature data for 19,549 proteins. In addition to focusing on the PTM sites themselves, residues spatially adjacent to the PTM sites in the AlphaFold2 structure were further calculated. The calculation method was to map the collected PTM sites into the AlphaFold2 structure and use Python to calculate the shortest spatial distance between the β-carbon atom (or α-carbon atom if it is glycine) of each non-PTM residue in the structure and the PTM site based on the residue atomic coordinates in the pdb file. Residues within a spatial distance of 4 Å (PTM_4Å) or 8 Å (PTM_8Å) from the PTM sites were also used as new PTM features.

[0092] Therefore, the sixth to eighth residue-level features are that the mutation occurs at a PTM site, the mutation occurs at a residue within a spatial distance of 4 Å from the PTM site, and the mutation occurs at a residue within a spatial distance of 8 Å from the PTM site, respectively.

[0093] Disordered regions: Experimentally verified protein disordered regions were obtained from two databases, MobiDB and DisProt, and their disordered regions were merged after removing redundancy, resulting in 2,530,679 disordered site feature data for 7754 proteins. Then, the ninth residue-level feature is that the mutation occurs at a disordered site.

[0094] Short linear motif: 16,687 short linear motif sites of 1,456 experimentally verified proteins were obtained from the Eukaryotic Linear Motif resource (ELM) database. The 10th residue-level feature was that mutations occurred at short linear motif sites.

[0095] Allosteric site: Experimentally confirmed PDB (Protein Data Bank) allosteric sites were obtained from the ASD database (Allosteric Database). Subsequently, these sites were checked for their match with the actual PDB structure sites. After removing the mismatched sites, these valid allosteric sites were mapped to the UniProt protein canonical sequences.

[0096] Since it was found that the collected experimental data was relatively limited, allosteric site data predicted by the AlloSitePro algorithm was also downloaded from ASD. The prediction was based on PDB and AlphaFold2 structures, and the two types of data were processed separately: (1) For allosteric sites predicted based on the PDB structure: Sites inconsistent with the actual PDB structure were removed, and the remaining allosteric sites were mapped to the UniProt canonical protein sequences; (2) For allosteric sites predicted based on the AlphaFold2 structure: The median pLDDT of all residues in the AlphaFold2 structure was calculated, and the allosteric sites predicted in the structure with a median pLDDT ≥ 70 were mapped to the UniProt canonical protein sequences. The allosteric site data from these three sources was integrated to obtain 1,435,333 allosteric site feature data of 14,160 proteins. The 11th residue-level feature was whether mutations occurred at allosteric sites.

[0097] Protein-protein binding sites: Experimental and predicted protein-protein binding site data were collected or calculated from the following two sources: (1) Experimental source data: This part of the data comes from protein homo / heterodimer crystal structures. Mutations were mapped to the corresponding structural chains, and while ensuring a high coverage of the mutant chains, as many dimer structures as possible were retained. Based on this, 5,623 homodimer structures of 3,166 proteins and 2,445 heterodimer structures of 1,708 proteins were obtained. For these dimer structures, the protein binding sites on the mutant chains were retrieved using the EBI-PISA API (a program interface name provided by EBI that provides protein interaction interfaces), and then the structures and protein standard sequences were mapped through sequence alignment, resulting in 164,891 protein-protein binding sites of 3,906 proteins; (2) Predicted data from Interactome INSIDER (a protein interaction annotation database): All human protein-protein interaction data were downloaded, and all interaction partner-specific interfacial residues of each protein were integrated as protein-protein binding sites, resulting in 1,221,720 protein-protein binding sites of 14,722 proteins, and the binding sites that could be mapped to high-quality AlphaFold2 structures were retained. Finally, the experimental and predicted sites were integrated to obtain 992,348 high-quality protein-protein binding site feature data of 10,667 proteins. Therefore, the 12th residue-level feature is that the mutation occurs at the protein-protein binding site.

[0098] Protein-small molecule ligand binding sites: Experimental and predicted protein-small molecule ligand binding site data were obtained from two sources: (1) Experimental data: This part of the data comes from the PDB crystal structure of the protein-small molecule ligand complex where the mutation occurs. Structures with a relatively high coverage of the mutant chain were also considered. A total of 35,404 mutant chains interacting with ligands were obtained, involving 2,934 proteins. The protein-ligand binding sites in the mutant chains were extracted through the EBI-PISA API, and sequence alignment was performed between the structural chain and the standard protein sequence to obtain 137,726 binding sites for 2,934 proteins. (2) P2Rank (name of the computational tool) predicted sites: P2Rank was used to predict potential protein-small molecule ligand binding sites in the AlphaFold2 structure. It identifies ligand-binding capabilities by analyzing the soluble local chemical neighborhood on the protein surface, obtaining 261,960 binding sites for 12,310 proteins. The experimental and predicted data were integrated to obtain 371,359 protein-small molecule ligand binding site feature data for 12,860 proteins. Therefore, the 13th residue-level feature is that the mutation occurs at the protein-small molecule ligand binding site.

[0099] Protein-nucleic acid binding sites: Protein-nucleic acid binding site data were collected from two sources: (1) Experimental data: From the PDB crystal structures of protein-nucleic acid ligand complexes where the mutations occurred. Structures with a relatively high coverage of the mutant chains were also considered. A total of 3,714 mutant chains interacting with nucleic acids were obtained, involving 394 proteins. Through the EBI-PISA API, the protein-nucleic acid binding sites in the mutant chains were extracted, and the structural sites were mapped to the sequences, resulting in 15,526 nucleic acid binding sites for 394 proteins. (2) GraphBind (a computational tool name) prediction data: 2,070 and 2,911 potential nucleic acid-binding proteins were collected from the ENPD (Eukaryotic nucleic acid binding protein database) and EuRBPDB (eukaryotic RNA binding proteins (RBPs) database) respectively. After merging and removing duplicates, the number of proteins was 4,245, and 2,521 of them had high-quality AlphaFold2 structures. GraphBind was used to predict the protein-nucleic acid binding sites of these AlphaFold2-structured proteins, obtaining 251,885 potential nucleic acid binding sites for 2,511 proteins. However, according to the subsequent results, only the experimental data were retained as the final protein-nucleic acid binding site characteristic data. Therefore, the 14th residue-level feature is that the mutation occurs at the protein-nucleic acid binding site.

[0100] Protein domains: Protein domain data were obtained from InterPro (v101.0) (InterPro: A Database of Protein Families, Domains and Functional Sites). There were 4,396,188 domain sites for 13,505 proteins. Then the 15th residue-level feature is that the mutation occurs at the protein domain site.

[0101] Physicochemical properties of amino acid residues (6 groups): According to the basic physicochemical properties of 20 standard amino acids, they were divided into 6 groups: Aliphatic amino acids include alanine, isoleucine, leucine, methionine, and valine; Aromatic amino acids include phenylalanine, tryptophan, and tyrosine; Positively charged amino acids are histidine, lysine, and arginine; Negatively charged amino acids include aspartic acid and glutamic acid; Neutral amino acids include asparagine, glutamine, serine, and threonine; Special amino acids include cysteine, proline, and glycine.

[0102] Therefore, the 16th to 21st residue-level features are that the mutant reference residues belong to aliphatic (alanine, isoleucine, leucine, methionine, and valine), the mutant reference residues belong to aromatic (phenylalanine, tryptophan, and tyrosine), the mutant reference residues belong to positively charged (histidine, lysine, and arginine), the mutant reference residues belong to negatively charged (aspartic acid and glutamic acid), the mutant reference residues belong to neutral (asparagine, glutamine, serine, and threonine), and the mutant reference residues belong to special (cysteine, proline, and glycine).

[0103] In this embodiment, a total of 6 mutant-level features are defined and calculated as follows:

[0104] Change in the charge of amino acid residues before and after mutation (3 groups): Charge values are assigned to 20 standard amino acid residues, where lysine and arginine are defined as positively charged (+1), aspartic acid and glutamic acid are defined as negatively charged (-1), and other amino acids are defined as neutral (0). According to the absolute value of the change in the charge of amino acid residues before and after mutation, the mutations are divided into three categories: Charge_0 (charge change is 0), Charge_1 (charge change is 1), and Charge_2 (charge change is 2).

[0105] Therefore, the 1st to 3rd mutant-level features are that the absolute value of the change in the charge of amino acid residues before and after mutation is 0, the absolute value of the change in the charge of amino acid residues before and after mutation is 1, and the absolute value of the change in the charge of amino acid residues before and after mutation is 2.

[0106] Change in the volume of amino acid residues before and after mutation (3 groups): The 20 standard amino acid residues are divided into three groups according to their volume size, namely: small group (alanine, glycine, serine, cysteine, proline, threonine, aspartic acid, and asparagine), encoded as 0; medium group ( histidine, glutamic acid, and glutamine), encoded as 1; large group (isoleucine, leucine, methionine, lysine, arginine, phenylalanine, tryptophan, and tyrosine), encoded as 2. According to the absolute value of the change in the volume of amino acid residues before and after mutation, the mutations are divided into three categories: Size_0 (volume change is 0), Size_1 (volume change is 1), and Size_2 (volume change is 2).

[0107] Therefore, the 4th to 6th mutant-level features are that the absolute value of the change in the volume of amino acid residues before and after mutation is 0, the absolute value of the change in the volume of amino acid residues before and after mutation is 1, and the absolute value of the change in the volume of amino acid residues before and after mutation is 2.

[0108] For the above 21 residue-level features and 6 mutation-level features, the evaluation values of each feature are calculated using statistical metrics. In this embodiment, the statistical metric used is the odds ratio. The odds ratio (Odds Ratio, OR) is a statistical metric used to evaluate the likelihood of an event occurring in one group relative to another group. The association differences between pathogenic mutations and benign mutations in the 21 residue-level features and 6 mutation-level features are quantified through the odds ratio. The data calculation expression for OR is:

[0109]

[0110] Among them, represents the odds ratio of pathogenic missense mutations and benign missense mutations in a molecular phenotypic feature, represents the number of pathogenic missense mutations that conform to this molecular phenotypic feature, represents the number of pathogenic missense mutations that do not conform to this molecular phenotypic feature, represents the number of benign missense mutations that conform to this molecular phenotypic feature, represents the number of benign missense mutations that do not conform to this molecular phenotypic feature.

[0111] The 27 molecular phenotypic features include both features based on the standard sequence, such as PTM sites, disordered regions, short linear motifs, protein domains, six physicochemical properties, three changes in residue charge before and after mutation, and three changes in residue volume before and after mutation, and features based on the structure, such as core, surface, three secondary structures, PTM_4Å, PTM_8Å, allosteric sites, protein-protein binding sites, protein-small molecule ligand binding sites, and protein-nucleic acid binding sites.

[0112] When calculating, it should be noted that: (1) It is only carried out for proteins with current feature annotations. For example, when calculating the OR of PTM sites, only proteins with PTM sites in the sequence are considered; (2) For features based on the standard sequence, the analyzed range in each protein should be the entire standard sequence. (3) For structural features obtained from the AlphaFold2 predicted structure or crystal structure, the analysis range of each protein is limited to the standard sequence region that can be covered by the structure. For example, when performing the OR comparison analysis of PTM_4Å, proteins with AlphaFold2 structures in the dataset should be considered, and whether the standard sequence region covered by the AlphaFold2 structure has PTM_4Å residues should be concerned to determine whether the protein is included in the analysis.

[0113] In this embodiment, the two-sided Fisher's exact test is used to calculate the OR. This test is implemented through the fisher_exact() function in the Python SciPy library. The , , and These four statistics are used as inputs, and the odds ratio (OR) and p-value are calculated based on the two-sided Fisher's exact test method through the fisher_exact() function.

[0114] To control the influence of multiple hypothesis testing, the Benjamini-Hochberg (BH) multiple testing correction method is used to correct the p-value obtained from Fisher's exact test, and the corrected q-value is obtained. The 95% confidence interval is calculated using the Wald method, which provides an estimated range within which the true OR is likely to fall at the 95% confidence level. The confidence interval does not have a fixed numerical threshold, and its interpretation needs to consider the research background, sample size, effect size, etc. Generally, if the OR value is not equal to 1, it indicates that this feature can distinguish pathogenic and benign mutations to a certain extent. On this basis, if the 95% confidence interval of the OR value does not contain 1 and the q-value < 0.05, it indicates that this feature can significantly distinguish pathogenic and benign mutations. Among them, OR > 1 indicates that this feature is significantly associated with pathogenic mutations; OR < 1 indicates that this feature is significantly associated with benign mutations.

[0115] The results are as Figure 2 shown. The bar chart represents the magnitude of log 10 OR Pathogenic / Benign , the color represents different features, and the error bars represent the 95% confidence interval of log 10 OR Pathogenic / Benign . The smaller the range of the error bars, the smaller the range within which log 10 OR Pathogenic / Benign is likely to fall, and the more accurate and reliable the result. The number of "*" is used to show the p-value (q-value) of Fisher's exact test after BH correction, representing the statistical significance level of the result. Specifically: q-value < 0.05, "*"; q-value < 0.001, "**"; q-value < 0.0001, "***". If log 10 OR Pathogenic / Benign is not equal to 0, it indicates that this feature can distinguish pathogenic and benign mutations to a certain extent. On this basis, if log 10 OR Pathogenic / BenignThe 95% confidence interval does not include 0 (i.e., the 95% confidence interval of the OR described above does not include 1), and the q-value < 0.05, indicating that this feature can significantly distinguish pathogenic mutations from benign mutations. Charge_0, Charge_1, Charge_2 represent the changes in residue charge before and after mutation; Size_0, Size_1, Size_2 represent the changes in residue volume before and after mutation; PTM represents the residue at the post-translational modification site, PTM_4Å represents the residue at the post-translational modification site and the residues within a spatial distance of 4 Å therefrom, and PTM_8Å represents the residue at the post-translational modification site and the residues within a spatial distance of 8 Å therefrom.

[0116] As can be seen from the figure, among the 27 molecular phenotypic features, 24 features can significantly distinguish pathogenic mutations from benign mutations. Although the remaining 3 features are not statistically significant, they also show a certain degree of discrimination. Therefore, the log 10 OR Pathogenic / Benign of these 27 molecular phenotypic features is used as the subsequent modeling feature.

[0117] The mutation level features also include the changes in protein stability and protein-protein binding affinity caused by the mutation. Among them, the change in protein stability is used to evaluate the impact of the mutation on protein stability. In this example, the PremPS tool (name of the calculation tool) is used to calculate the change in structural stability caused by the mutation based on the AlphaFold2 structure. The calculated result (the change in free energy before and after the mutation calculated by PremPS, which is the change in protein stability, unit: kcal / mol) is used to represent the change in conformational stability free energy caused by a single missense mutation. When it represents that the mutation increases the structural stability, represents that the mutation reduces the structural stability.

[0118] The change in protein-protein binding affinity is used to evaluate the impact of the mutation on protein-protein binding affinity. In this example, the computational tool MutaBind2 is used to calculate the change in protein binding affinity caused by a single missense mutation for the heterodimer structure and its mutation. The calculation result is represented by (the change in free energy before and after the mutation calculated by MutaBind2, which is the change in protein binding affinity, unit: kcal / mol), indicating that the mutation increases the binding affinity between the interacting parties, while indicating that the binding affinity between the two is reduced, and the protein-protein interaction is disrupted to a certain extent.

[0119] Combined with Figure 5 , Figure 5Protein stability changes caused by pathogenic missense mutations (Pathogenic) and benign missense mutations (Benign) ( ), and changes in binding affinity ( ). Figure 5 In or To measure whether there is a significant difference between the two mutations, the p-value is obtained by t-test and the q-value is obtained by BH correction, and the significance level is divided: q-value ≥ 0.05, there is no significant difference between the two; q-value < 0.05, there is a significant difference between the two. Specifically: q-value < 0.05, "*"; q-value < 0.001, "**"; q-value < 0.0001, "***". As can be seen from Figure 5 , compared with benign mutations, pathogenic mutations can significantly cause a decrease in protein stability and protein-protein binding affinity, while the changes in stability and binding affinity caused by benign mutations are relatively small. Therefore, the protein stability changes and protein-protein binding affinity changes caused by mutations are also used as modeling features.

[0120] Next, the comprehensive evaluation values of the five major types of protein-level features are obtained by using the statistical index odds ratio OR for enrichment analysis and annotation of protein-level features. Specifically:

[0121] Protein-level features include 5 major categories, namely 4 major categories of protein feature function classifications at the system level, subcellular localization, KEGG (Kyoto Encyclopedia of Genes and Genomes) pathways, and NetSlim (a database name) cancer signaling pathways, as well as specific protein domains.

[0122] Protein function classification: 24 functional classifications of 13,987 human proteins were obtained from PANTHER (v19.0) (Protein ANalysis THrough Evolutionary Relationships), and this classification was obtained based on protein evolutionary relationships and ancestrally highly conserved functions. The 24 functional classifications are called 24 sub-features, and it is determined and analyzed whether the protein where the mutation is located conforms to the 24 sub-features respectively.

[0123] Subcellular localization: Human protein subcellular localization obtained from immunofluorescence staining experiments in the Human Protein Atlas (HPA) was used. Based on the reliability of experimental evidence from high to low, HPA classifies subcellular localization into four levels: Enhanced, Supported, Approved, and Uncertain. After excluding subcellular localizations at the Uncertain level, 36 subcellular localization information of 11,926 proteins was collected. In this example, the 36 subcellular localization information is called 36 subclass features, and it was respectively determined and analyzed whether the protein where the mutation is located conforms to these 36 subclass features.

[0124] KEGG pathways: 186 KEGG pathways and the genes in each pathway were obtained from the MSigDB database (Molecular Signatures Database), and all gene symbols were mapped to UniProt standard protein names. In this example, the 186 KEGG pathways are called 186 subclass features, and it was respectively determined and analyzed whether the protein where the mutation is located conforms to these 186 subclass features.

[0125] NetSlim cancer signaling pathways: 32 manually collected high-quality cancer signaling pathways and the genes in each pathway were obtained from the NetSlim database, and all gene symbols were mapped to UniProt standard protein names. In this example, the 32 NetSlim cancer signaling pathways are called 32 subclass features, and it was respectively determined and analyzed whether the protein where the mutation is located conforms to these 32 subclass features.

[0126] Protein specific domains: A total of 4,949 types of human protein specific domains were obtained from InterPro. In this example, the 4,949 types of specific domains are called 4,949 subclass features, and it was respectively determined and analyzed whether the protein where the mutation is located conforms to these 4,949 subclass features.

[0127] The calculation formula for the odds ratio of pathogenic missense mutations and benign mutations in subclass protein level features is:

[0128]

[0129] Among them, represents the odds ratio of pathogenic missense mutations and benign missense mutations in a subclass protein level feature, represents the number of pathogenic missense mutations that conform to this subclass protein level feature, Indicates the number of pathogenic missense mutations that do not conform to the protein level characteristics of this subclass, Indicates the number of benign missense mutations that conform to the protein level characteristics of this subclass, Indicates the number of benign missense mutations that do not conform to the protein level characteristics of this subclass.

[0130] For the protein characteristics at the four major system levels, the background includes the complete sequences of all proteins in the pathogenic mutation dataset or the benign mutation dataset; for a specific domain, the background includes the complete sequences of all proteins containing this domain in the pathogenic mutation dataset or the benign mutation dataset.

[0131] In summary, the OR of each subclass characteristic can be calculated. pathogenic / benign Then, for any mutation, the formula for calculating the comprehensive value of the odds ratio of each of the five major protein level characteristics is as follows:

[0132]

[0133] Where, Indicates any one of the five major protein level characteristics: protein function classification, human subcellular localization, KEGG pathway, NetSlim cancer signaling pathway, and specific protein domain, Indicates the comprehensive value of the odds ratio of pathogenic missense mutations and benign missense mutations in one major protein level characteristic, Is the number of all subclass characteristics in the major protein characteristic that the protein where the mutation is located conforms to, Is the Odds ratio of the th subclass characteristic.

[0134] The missense mutations in the training set and the OR evaluation values of their corresponding 27 molecular phenotype characteristics, the comprehensive OR evaluation values of each of the 5 major protein level characteristics, and the predicted values of the changes in protein stability and binding affinity caused by the mutation are used to form the feature matrix of the model training dataset.

[0135] Finally, a random forest classification model is constructed. Using the feature matrix of the model training dataset as input, the optimal hyperparameter combination is determined through the system's hyperparameter optimization strategy to maximize the prediction performance of the model. Specifically, two key hyperparameters are tuned: n_estimators (the number of decision trees in the random forest) and max_features (the maximum number of features in the random forest). Among them, n_estimators represents the number of decision trees in the random forest, with a value range from 50 to 1950 and a step size of 50. This range covers forest configurations from small to large scale to comprehensively evaluate the impact of different numbers of trees on the model performance. And max_features controls the maximum number of features considered by each tree during construction, and its values include strategies of using all features (None), selecting features by the square root of the features (sqrt), and selecting features by the logarithm (log2). To comprehensively explore the configuration space of these hyperparameters, a parameter grid containing all possible parameter combinations is constructed, and on this basis, five-fold cross-validation of the training set is carried out. This strategy ensures that the impact of each combination on the model performance can be systematically evaluated, so as to select the optimal hyperparameter combination. Finally, when the number of decision trees is equal to 1950 and the strategy of selecting features by logarithm is used, the best average Matthews Correlation Coefficient (MCC) of five-fold cross-validation is obtained, and the missense mutation functional effect prediction model MutaPheno is obtained.

[0136] Next, the performance of the missense mutation functional effect prediction model and other existing models will be evaluated and compared using an independent test set.

[0137] To compare the prediction effects of different models, the following evaluation metrics are introduced: Recall, that is, the True Positive Rate (TPR), which is used to represent the proportion of positive mutations that the model can correctly identify, that is, the proportion of positive mutations predicted as positive among positive mutations. Its expression is: ;

[0138] Among them, represents Recall, represents the true positives, that is, the number of mutations that are actually positive and predicted as positive, represents the false negatives, that is, the number of mutations that are actually positive and predicted as negative.

[0139] The False Positive Rate (FPR) represents the proportion of negative mutations misjudged as positive mutations. Its expression is: ;

[0140] Among them, represents the false positive rate, represents false positives, that is, the number of mutations that are actually negative but predicted to be positive, represents true negatives, that is, the number of mutations that are actually negative and predicted to be negative.

[0141] Precision, which represents the proportion of true positive mutations among the mutations predicted to be positive by the model, that is, the proportion of positive mutations among the mutations predicted to be positive by the model, and its expression is: ;

[0142] where, represents precision, represents true positives, that is, the number of mutations that are actually positive and predicted to be positive, represents false positives, that is, the number of mutations that are actually negative but predicted to be positive.

[0143] F1-score, which comprehensively considers the precision and recall metrics of the model, is the harmonic mean of the two, and is applicable to scenarios with imbalanced positive and negative data, and its expression is: ;

[0144] where, represents the F1 score value, represents recall, represents precision.

[0145] MCC (Matthews Correlation Coefficient), which is also an indicator for comprehensively evaluating the performance of the model and is reliable in the case of imbalanced positive and negative data, and its expression is: ;

[0146] where, represents true positives, that is, the number of mutations that are actually positive and predicted to be positive, represents true negatives, that is, the number of mutations that are actually negative and predicted to be negative, represents false positives, that is, the number of mutations that are actually negative but predicted to be positive, represents false negatives, that is, the number of mutations that are actually positive and predicted to be negative.

[0147] AUROC (Area Under the Receiver Operating Characteristic Curve), which is a commonly used metric for evaluating the performance of classification models, represents the ability of the model to distinguish between positive and negative samples. AUROC is the area under the ROC curve (Receiver Operating Characteristic), and the closer the value is to 1, the stronger the classification ability of the model and the better it can distinguish between positive and negative mutations.

[0148] In this embodiment, an independent test set of pathogenic mutations was first constructed: pathogenic and benign missense mutations in the AlphaMissense and CPT-1 test sets were downloaded, both of which were from the ClinVar database. The identifiers of the mutations in the AlphaMissense test set are UniProt standard proteins, and the identifiers provided for the mutations in the CPT-1 test set are Gene Symbols. Mutation standardization and protein structure mapping were performed on the mutation data from both sources, and the operations were the same as those described above, so no further elaboration will be given here. The standardized mutation data were merged, only the missense mutations successfully mapped to the structure were retained, and the mutations repeated with the training set were removed to ensure the independence of the test set. Finally, 28,717 protein missense mutations were obtained as an independent test set of pathogenic mutations, including 12,539 and 16,178 pathogenic / benign mutations. In this test, 49 mutation functional effect prediction models were used, including 47 existing models in the dbNSFP database (Database for Non-Synonymous Functional Polymorphisms), the deep learning model CPT-1, and the missense mutation functional effect prediction model MutaPheno of the present invention. The above 49 models were tested using the independent test data set of pathogenic mutations. The initial evaluation results showed that the performance of some models was extremely superior. After in-depth exploration, it was found that there were mainly two reasons: one was data leakage - 12 models (such as MetaRNN, ClinPred, etc.) used mutations in ClinVar during training, and these mutations were also included in the test set, resulting in the artificial elevation of the model performance; the other was that some models only provided prediction results for a very small number of mutations in the test set and could not comprehensively reflect their overall performance level (for example, two CADD models only made predictions for 565 mutations). Therefore, to ensure the fairness and objectivity of the evaluation, these 14 models were excluded from this comparative analysis.

[0149] Among the remaining 35 models, the model MutaPheno proposed by the present invention ranked sixth in terms of the AUROC metric, demonstrating competitive and robust performance. On the complete test set (n = 28,717), the top ten models were: AlphaMissense, CPT-1, VEST4, DEOGEN2, MetaSVM, MutaPheno, M_CAP, MVP, Eigen, and MutFormer.

[0150] To further ensure the fairness of the evaluation and eliminate the bias caused by differences in mutation predictability among different tools, the mutation intersections of the prediction results of the top 10 models were calculated. Finally, a subset containing 23,222 mutations was obtained, including 10,485 pathogenic mutations and 12,737 benign mutations. On this shared mutation set, the performance of each model was evaluated based on multiple evaluation metrics, including AUROC, Matthews correlation coefficient (MCC), F1-score, Precision, and Recall, and the test results are shown in Table 1.

[0151] Table 1: Performance of the top 10 models on the independent test set of pathogenic mutations

[0152] Notably, on this shared test set, MutaPheno ranked seventh overall, achieving an area under the receiver operating characteristic curve (AUROC) of 0.883, an F1-score of 0.794, and an MCC of 0.619, demonstrating balanced and excellent prediction performance for both types of mutations. Figure 3 The ROC curves of the top 7 models in terms of AUROC on the pathogenic shared test set are shown, and the ROC curve plots the relationship between TPR and FPR.

[0153] To further evaluate the robustness and generalization ability of MutaPheno in identifying cancer driver mutations, an independent cancer driver mutation test data set was constructed, and the performance of the model was evaluated on this test set.

[0154] For the cancer driver mutation test set, driver mutation and passenger mutation data were collected from three different sources, namely: the OncoKB database (Oncology Knowledge Base), the CGI (Cancer Genome Interpreter) platform, and the comprehensive cancer missense mutation data set (MutaGene).

[0155] OncoKB is an evidence-based precision oncology knowledge base that provides researchers with detailed somatic mutation annotation information, including the biological effects and cancer driver status of mutations. Missense mutations annotated as "oncogenic / likely oncogenic" or "likely neutral" are retained. The identifiers provided by OncoKB mutations are genomic database (Ensembl) transcripts.

[0156] CGI is a multifunctional platform for assessing and grading the driver effects of tumor mutations based on existing knowledge and evidence. It integrates data from multiple databases and literature sources. The identifiers provided by CGI mutations are Ensembl transcripts.

[0157] The comprehensive cancer missense mutation dataset is integrated from multiple sources and annotated based on experimental validation results of the characteristics of mutations such as protein function, binding, and transformation ability. Thus, the cancer missense mutations therein can be classified into two categories: "non-neutral" or "neutral". The identifiers provided by the mutations are gene symbols (Gene Symbol).

[0158] The "likely oncogenic" and "oncogenic" mutations in OncoKB, all oncogenic missense mutations in CGI, and the "non-neutral" missense mutations in the comprehensive dataset are classified as driver mutations; while the "likely neutral" mutations in OncoKB and the "neutral" missense mutations in the comprehensive cancer missense mutation dataset (MutaGene) are classified as passenger mutations.

[0159] Mutation standardization and protein structure mapping are performed on the mutation data from the three sources, and the operations are the same as above, so they will not be elaborated here. After merging the standardized mutations, only the missense mutations that are successfully mapped to the structure are retained, and the mutations that are repeated with the training set are removed to ensure the independence of the test set. Finally, 6,322 missense mutations are obtained as an independent test set for driver mutations, including 2,317 driver mutations and 4,005 passenger mutations.

[0160] The performance of the model was evaluated on the above independent driver mutation test dataset. In this benchmark test, a total of 53 models were compared, covering 47 prediction tools from dbNSFP, CPT-1, MutaPheno of this embodiment, and 4 tools specifically for predicting cancer driver mutations: CHASM, CHASMplus, ParsSNP, and AI-Driver.

[0161] On the complete test set, MutaPheno performed best overall, outperforming all other methods in terms of the AUROC metric and showing strong performance on multiple evaluation metrics.

[0162] For a more rigorous and fair comparison, 12 representative models were selected for analysis, including: 6 models that performed better than MutaPheno in the pathogenic mutation sharing test set (based on the AUROC metric); the MVP model, which showed the best performance in predicting driver mutations in the dbNSFP tool; 4 tools specifically for cancer driver mutation prediction (CHASM, CHASMplus, ParsSNP, AI-Driver); and MutaPheno proposed by the present invention.

[0163] The mutation intersection with prediction results for all 12 models was extracted to obtain a shared dataset containing 5,770 missense mutations (2,018 of which were driver mutations and 3,752 were passenger mutations). The comparative performance of these 12 models is summarized in Table 2.

[0164] Table 2: Performance of 12 models on the independent test set of cancer driver mutations

[0165] MutaPheno achieved the highest area under the receiver operating characteristic curve (AUROC) on this dataset, which was 0.841, and maintained the overall first place, highlighting its excellent performance in distinguishing driver mutations from passenger mutations. Figure 4 The ROC curves of the top 7 models in terms of AUROC on the driver sharing test set are shown. The ROC curve plots the relationship between TPR and FPR.

[0166] In summary, the missense mutation functional effect prediction model MutaPheno proposed by the present invention, as an interpretable feature-based model, not only showed better performance in pathogenic mutation prediction but also performed best in the identification of driver mutations, highlighting the wide applicability of MutaPheno beyond general functional effect prediction. Although MutaPheno was not specifically trained on cancer-related datasets, it showed strong discrimination ability in identifying cancer driver mutations and even outperformed specialized cancer driver mutation prediction tools in some evaluation metrics. In contrast, other models such as AlphaMissense, CPT-1, and VEST4, although having certain advantages in pathogenic mutation prediction, performed poorly in distinguishing driver mutations. This finding indicates that the molecular and functional features integrated into the model have high information value in various pathogenic backgrounds including tumorigenesis.

[0167] Example 2:

[0168] Based on Example 1, this example also provides a method for predicting the functional effect of missense mutations based on molecular phenotype features, including:

[0169] Obtain the position, reference residue, and mutant residue of the missense mutation to be tested in the UniProt standard sequence; determine whether the missense mutation to be tested conforms to 27 molecular phenotype characteristics. If it does not conform to a certain molecular phenotype characteristic, the annotation result of this molecular phenotype characteristic is 0. If it conforms to a certain molecular phenotype characteristic, the annotation result of this molecular phenotype characteristic is the evaluation value of this molecular phenotype characteristic, and obtain the molecular phenotype characteristic judgment result.

[0170] Determine whether the missense mutation to be tested conforms to five major categories of protein-level characteristics. If the protein where the missense mutation to be tested is located does not conform to any subclass in a certain major category of protein-level characteristics, the annotation result of this protein-level characteristic is 0. If the protein where the missense mutation to be tested is located conforms to one or more subclasses in a certain major category of protein-level characteristics, the annotation result of this protein-level characteristic is the comprehensive evaluation value obtained by adding the evaluation values of these one or more subclasses, and obtain the five major categories of protein-level characteristic judgment results.

[0171] Use the PremPS and MutaBind2 tools to calculate the predicted values of the protein stability change and binding affinity change of the missense mutation to be tested.

[0172] Take the molecular phenotype characteristic judgment result, protein-level characteristic judgment result, and the predicted values of protein stability change and binding affinity change, a total of 34 annotation characteristics, as inputs, and based on the missense mutation functional effect prediction model MutaPheno described in Example 1, predict the functional effect of the missense mutation to be tested, and output the functional effect classification label of the missense mutation to be tested:

[0173] When the classification label is 1, it indicates that the missense mutation to be predicted is a mutation with functional impact. In the training stage, it corresponds to a pathogenic missense mutation, and in practical applications, it can be generalized to identify cancer driver mutations. When the classification label is 0, it indicates that the predicted mutation is a mutation with no obvious functional impact. In the training stage, it corresponds to a benign missense mutation, and in practical applications, it can be generalized to identify cancer passenger mutations.

[0174] The above is only the preferred embodiment of the present invention. It should be noted that for those of ordinary skill in the art, without departing from the technical principle of the present invention, several improvements and deformations can be made, and these improvements and deformations should also be regarded as the protection scope of the present invention.

Claims

1. A method for constructing a prediction model for the functional effect of a missense mutation, characterized in that: include: Integrating missense mutation data from different sources to construct a training dataset, wherein the missense mutations include pathogenic missense mutations and benign missense mutations; For each missense mutation in the training data set, molecular phenotypic characteristics and protein-level characteristics were enriched and analyzed to obtain the molecular phenotypic characteristics evaluation value and the comprehensive evaluation value of the major protein-level characteristics, and the predicted value of the protein stability change and binding affinity change caused by the missense mutation was calculated. According to the molecular phenotype feature evaluation value and the comprehensive evaluation value of the large-category protein level feature, each missense mutation in the training data set is annotated with a molecular phenotype feature and annotated with a large-category protein level feature to obtain a molecular phenotype feature annotation result and a large-category protein level feature annotation result, respectively; A multidimensional feature vector is formed using the molecular phenotype feature annotation results of each missense mutation in the training data set, the feature annotation results of the major protein level, and the predicted values ​​of the protein stability change and binding affinity change caused by the missense mutation; Build a random forest classification model; The random forest classification model was trained with the multidimensional feature vector as input to obtain a missense mutation functional effect prediction model.

2. The method for constructing a missense mutation functional effect prediction model according to claim 1, characterized in that: The missense mutation data includes an identifier, a missense mutation position, a reference residue, and a mutated residue.

3. The method for constructing a missense mutation functional effect prediction model according to claim 1, characterized in that: After constructing the training dataset, it also includes missense mutation normalization and protein structure mapping of the training dataset; The missense mutation standardization is used to unify missense mutation data using the standard protein name of the unified protein resource database as an identifier, specifically including: When the identifier of the missense mutation is a gene symbol, the gene symbol where the missense mutation is located is converted into a standard protein name of the unified protein resource database, and the position of the missense mutation and the reference residue are checked to see whether they are consistent with the residue at the same position of the standard protein sequence in the unified protein resource database. If they are consistent, the missense mutation is retained; if not, the missense mutation is deleted; When the identifier of the missense mutation is a standard protein name in the unified protein resource database, the missense mutation position and the reference residue are checked to see if they are consistent with the residue at the same position of the standard protein sequence in the unified protein resource database. If they are consistent, the missense mutation is retained; if they are inconsistent, the missense mutation is deleted. When the identifier of the missense mutation is a transcript in the reference sequence database, a transcript in the genome database, or a protein in the genome database, the transcript in the reference sequence database, the transcript in the genome database, or the protein in the genome database where the missense mutation is located is converted into the name of the isoform protein in the unified protein resource database, and the position of the missense mutation and the reference residue are checked to see whether they are consistent with the residue at the same position of the isoform protein sequence in the unified protein resource database. If they are consistent, the mutation is retained; if not, the mutation is deleted; Perform a double sequence alignment of the standard protein sequence of the unified protein resource database and the isoform protein sequence of the unified protein resource database, map the mutation position of the isoform protein sequence originally based on the unified protein resource database to the mutation position of the standard protein sequence based on the unified protein resource database, and check whether the reference residue of the missense mutation is consistent with the residue of the corresponding missense mutation position of the standard protein sequence of the unified protein resource database, if consistent, retain it, if inconsistent, delete the missense mutation; The retained missense mutations were used as missense mutations after normalization; The protein structure mapping comprises: According to the position and reference residues of the standardized missense mutation in the standard protein sequence of the unified protein resource database, the missense mutation is mapped to the AlphaFold2 structural sequence. If the standard protein sequence of the unified protein resource database of the protein where the missense mutation is located is consistent with the AlphaFold2 structural sequence, it is directly mapped. If the standard protein sequence of the unified protein resource database of the protein where the missense mutation is located is inconsistent with the AlphaFold2 structural sequence, a double sequence comparison is performed between the standard protein sequence of the unified protein resource database and the AlphaFold2 structural sequence to complete the mapping.

4. The method for constructing a missense mutation functional effect prediction model according to claim 1, characterized in that: The enrichment analysis of the molecular phenotype characteristics is performed on each missense mutation data in the training data set to obtain the molecular phenotype characteristic evaluation value, including: Calculate the number of pathogenic missense mutations and benign missense mutations that meet the molecular phenotype characteristics and the number of pathogenic missense mutations and benign missense mutations that do not meet the molecular phenotype characteristics in the training data set; The number of pathogenic missense mutations and benign missense mutations that meet the molecular phenotypic characteristics and the number of pathogenic missense mutations and benign missense mutations that do not meet the molecular phenotypic characteristics are used as input, and the odds ratio of pathogenic mutations and benign mutations in the molecular phenotypic characteristics is calculated based on the two-sided Fisher test method; The logarithmic function of the odds ratio of the molecular phenotypic feature with base 10 was used as the evaluation value of the molecular phenotypic feature.

5. The method for constructing a missense mutation functional effect prediction model according to claim 4, characterized in that: The molecular phenotypic characteristics are divided into residue level characteristics and mutation level characteristics; The residue level features include that the residue where the mutation occurs is a core residue, the residue where the mutation occurs is a surface residue, the secondary structure of the protein where the mutation occurs is a helical type, the secondary structure of the protein where the mutation occurs is a folding type, the secondary structure of the protein where the mutation occurs is a ring structure type, the mutation occurs at a post-translational modification site, the mutation occurs at a residue within a spatial distance of 4Å from the post-translational modification site, the mutation occurs at a residue within a spatial distance of 8Å from the post-translational modification site, the mutation occurs at a disordered region site, the mutation occurs at a short linear motif site, the mutation occurs at an allosteric site, the mutation occurs at a protein-protein binding site, the mutation occurs at a protein-small molecule ligand binding site, the mutation occurs at a protein-nucleic acid binding site, the mutation occurs at a protein domain site, the physicochemical properties of the mutation reference residue are aliphatic, the physicochemical properties of the mutation reference residue are aromatic, the physicochemical properties of the mutation reference residue are positively charged, the physicochemical properties of the mutation reference residue are negatively charged, the physicochemical properties of the mutation reference residue are neutral, and the physicochemical properties of the mutation reference residue are special; The mutation level characteristics include the absolute value of the change in charge of the amino acid residues before and after the mutation being 0, the absolute value of the change in charge of the amino acid residues before and after the mutation being 1, the absolute value of the change in charge of the amino acid residues before and after the mutation being 2, the absolute value of the change in volume of the amino acid residues before and after the mutation being 0, the absolute value of the change in volume of the amino acid residues before and after the mutation being 1, and the absolute value of the change in volume of the amino acid residues before and after the mutation being 2.

6. The method for constructing a missense mutation functional effect prediction model according to claim 5, characterized in that: The physicochemical properties of the mutation reference residues are based on the classification of 20 standard amino acids, specifically including: Alanine, isoleucine, leucine, methionine and valine belong to the aliphatic group; phenylalanine, tryptophan and tyrosine belong to the aromatic group; histidine, lysine and arginine are positively charged; aspartic acid and glutamic acid are negatively charged; asparagine, glutamine, serine and threonine are neutral; cysteine, proline and glycine are special; The method for determining the absolute value of the charge change of the amino acid residues before and after the mutation includes: Set the charge parameters, including: lysine and arginine are +1 positive charge, aspartic acid and glutamic acid are -1 negative charge, and the remaining amino acids are 0 neutral charge; According to the charge parameter and the amino acid types of the amino acid residues before and after the mutation, the absolute value of the charge change of the amino acid residues before and after the mutation is calculated; The method for determining the absolute value of the volume change of the amino acid residue before and after the mutation includes: Set volume coding, specifically: the volume coding of alanine, glycine, serine, cysteine, proline, threonine, aspartic acid and asparagine is 0, the volume coding of valine, histidine, glutamic acid and glutamine is 1, and the volume coding of isoleucine, leucine, methionine, lysine, arginine, phenylalanine, tryptophan and tyrosine is 2; According to the volume code and the amino acid types of the amino acid residues before and after the mutation, the absolute value of the volume change of the amino acid residues before and after the mutation was calculated.

7. The method for constructing a missense mutation functional effect prediction model according to claim 1, characterized in that: The major protein-level features include human protein function categories, human protein subcellular localization information categories, standard protein ranges of the unified protein resource database corresponding to Kyoto Encyclopedia of Genes and Genomes pathways and cancer signaling pathways, and protein specific domain ranges; Each of the major protein-level features includes multiple minor protein-level features, and the comprehensive evaluation value of the major protein-level features is obtained by summing the evaluation values ​​of all minor protein-level features that are consistent with the protein where the missense mutation is located in each major protein-level feature.

8. The method for constructing a missense mutation functional effect prediction model according to claim 7, characterized in that: The enrichment analysis of protein-level features is performed on each missense mutation in the training data set to obtain a comprehensive evaluation value of a large category of protein-level features, including: The number of pathogenic missense mutations and benign missense mutations that meet the small class protein level characteristics and the number of pathogenic missense mutations and benign missense mutations that do not meet the small class protein level characteristics in the training data set are calculated, and used as input to calculate the odds ratio of each small class protein level characteristic based on the two-sided Fisher test method; For any missense mutation, the odds ratios of all small-category protein-level features that the protein containing the missense mutation meets in each large-category protein-level feature are summed using a logarithmic function with base 10, and this sum is used as the comprehensive evaluation value of the missense mutation in the corresponding large-category protein-level feature.

9. The method for constructing a missense mutation functional effect prediction model according to claim 4 or 8, characterized in that: The calculation expression of the odds ratio of missense mutations in molecular phenotypic characteristics is: ； in, represents the odds ratio of pathogenic missense mutations to benign missense mutations in a molecular phenotypic feature, Indicates the number of pathogenic missense mutations that meet the molecular phenotype characteristics. Indicates the number of pathogenic missense mutations that do not meet the molecular phenotype characteristics. Indicates the number of benign missense mutations that meet the molecular phenotype characteristics. represents the number of benign missense mutations that do not conform to the molecular phenotype characteristics; The calculation expression for the odds ratio of missense mutations in small class protein level features is: ； in, represents the odds ratio of pathogenic missense mutations to benign missense mutations in a small class of protein-level features, represents the number of pathogenic missense mutations that meet the protein level characteristics of this subclass, represents the number of pathogenic missense mutations that do not meet the protein level characteristics of this subclass, represents the number of benign missense mutations that meet the protein level characteristics of this subclass, represents the number of benign missense mutations that do not conform to the protein level characteristics of this subclass; The calculation expression of the comprehensive evaluation value of the protein level features of major categories is: ； in, represents any one of the broad classes of protein-level features, It represents the comprehensive value of the odds ratio of pathogenic missense mutations and benign missense mutations in this large category of protein-level features. is the number of all sub-category features that the mutated protein meets in the feature of this large category of proteins, It is one of the The odds ratio of a small class feature.

10. The method for constructing a missense mutation functional effect prediction model according to claim 1, characterized in that: According to the molecular phenotype feature evaluation value and the comprehensive evaluation value of the large-category protein level feature, each missense mutation in the training data set is annotated with a molecular phenotype feature and annotated with a large-category protein level feature to obtain a molecular phenotype feature annotation result and a large-category protein level feature annotation result, respectively, specifically including: For any missense mutation in the training data set, the molecular phenotypic features that match the missense mutation are annotated as the evaluation value of the corresponding molecular phenotypic features, and the molecular phenotypic features that do not match are annotated as 0, to obtain the molecular phenotypic annotation results; the subclass protein level features that match the missense mutation are annotated as the evaluation value of the corresponding subclass protein level features, and the subclass protein level features that do not match are annotated as 0, to obtain the annotation results of the subclass protein level features, and the annotation results of the subclass protein level features are summed up to obtain the annotation results of the large class protein level features.

11. A method for predicting the functional effect of a missense mutation, characterized in that: include: Obtaining the position, reference residue and mutant residue of the missense mutation to be tested in the standard sequence of the unified protein resource database; According to the position of the missense mutation to be tested in the standard sequence of the unified protein resource database, the reference residue and the mutant residue, it is judged whether the missense mutation to be tested conforms to the molecular phenotypic characteristics. If it does not conform to a certain molecular phenotypic characteristic, the annotation result of the molecular phenotypic characteristic is 0. If it conforms to a certain molecular phenotypic characteristic, the annotation result of the molecular phenotypic characteristic is the evaluation value of the molecular phenotypic characteristic, and the molecular phenotypic characteristic judgment result is obtained; According to the position, reference residue and mutant residue of the missense mutation to be tested in the standard sequence of the unified protein resource database, for each major category of protein level features, the evaluation values ​​of all the minor category of protein level features that the protein where the missense mutation to be tested is located in each major category of protein level features are summed to obtain a comprehensive evaluation value of the major category of protein level features, and the comprehensive evaluation value of all the major category of protein level features is used as the protein level feature judgment result; Calculate the predicted value of protein stability change and binding affinity change caused by the missense mutation to be tested; The molecular phenotype feature judgment results, protein level feature judgment results, and predicted values ​​of protein stability changes and binding affinity changes are used as input, and the missense mutation functional effect prediction model obtained based on the missense mutation functional effect prediction model construction method according to any one of claims 1 to 10 is used to predict the functional effect of the missense mutation to be tested, and the functional effect classification result of the missense mutation to be tested is obtained as output.

12. The method for predicting the functional effect of a missense mutation according to claim 11, characterized in that: The functional effect classification result includes: when the classification label is 1, it indicates that the missense mutation to be tested is a missense mutation with functional impact; when the classification label is 0, it indicates that the missense mutation to be tested is a missense mutation without functional impact; Among them, the missense mutation with functional impact includes pathogenic missense mutation or cancer-driving missense mutation; the missense mutation without functional impact includes benign missense mutation or passenger missense mutation.

Images