Method for constructing missense mutation functional effect prediction model and prediction method

Abstract

The present invention discloses a method for constructing a missense mutation functional effect prediction model and a prediction method, which relates to the field of bioinformatics technology. Based on pathogenic and benign missense mutation data, the present invention performs enrichment annotation on various molecular phenotypic features through statistical indicators to construct a machine learning model for predicting mutation functional effects. The model integrates feature information at multiple molecular levels, including comprehensive evaluation indicators at the residue level, mutation level, and protein level, and is trained and optimized based on the random forest algorithm. The present invention can not only accurately predict pathogenic mutations, but also exhibits excellent performance in the identification of cancer driver missense mutations, has 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.

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 the high cost and long cycle, it is difficult to be applied on a large scale for the screening of all potential mutations.

[0003] For this reason, various computational prediction methods and tools have been developed in recent years to evaluate the functional 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 the fitness-based scoring system and lack biological interpretability, making it difficult to reveal the underlying 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 is to generate a selective advantage through the acquisition of functional changes during the evolution of cancer cells. Therefore, the functional effect is the common basis of pathogenicity and drivability, and there is a significant overlap between these 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 the 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 strong 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 purpose, 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, enrichment analysis of molecular phenotypic characteristics and protein-level characteristics is performed to obtain the molecular phenotypic characteristic evaluation value and the comprehensive evaluation value of the major protein-level characteristics respectively, and the predicted values of the protein stability change and binding affinity change caused by the missense mutation are calculated;

[0010] According to the molecular phenotypic characteristic evaluation value and the comprehensive evaluation value of the major protein-level characteristics, molecular phenotypic characteristic annotation and major protein characteristic annotation are performed on each missense mutation in the training dataset to obtain the molecular phenotypic characteristic annotation result and the major protein-level characteristic annotation result respectively;

[0011] Using the molecular phenotypic characteristic annotation result, major protein-level characteristic 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] Using the multi-dimensional feature vector as the input, 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 of 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 to 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 they are inconsistent, 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 they are inconsistent, 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 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 molecular phenotype characteristics is performed on each missense mutation data in the training dataset to obtain the molecular phenotype characteristic evaluation value, including:

[0025] 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 dataset.

[0026] Take 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 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 logarithm function to the base 10 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 encoding, specifically including: the volume encoding of alanine, glycine, serine, cysteine, proline, threonine, aspartic acid and asparagine is 0, the volume encoding of valine, histidine, glutamic acid and glutamine is 1, and the volume encoding 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 encoding and the amino acid type of the amino acid residue before and after mutation.

[0039] Furthermore, the above-mentioned 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 pathway and the unified protein resource database corresponding to the cancer signaling pathway, and the protein specific domain range;

[0040] Each of the above-mentioned large-class protein level features 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 the small-class protein level features that the missense mutation in the protein conforms to in each large-class protein level feature.

[0041] Furthermore, 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 feature includes:

[0042] Calculate the number of pathogenic missense mutations and benign missense mutations that conform to the small-class protein level feature and the number of pathogenic missense mutations and benign missense mutations that do not conform to the small-class protein level feature 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 the 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] Furthermore, 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 one 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;

[0047] The calculation expression for the odds ratio of missense mutations in the feature at the subclass protein level is:

[0048] ;

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

[0050] The calculation expression for the comprehensive evaluation value of the feature at the superclass protein level is:

[0051] ;

[0052] Among them, represents any one of the features at the superclass protein level, represents the comprehensive value of the odds ratio of pathogenic missense mutations and benign missense mutations in this feature at the superclass protein level, is the number of all subclass features that the protein where the mutation is located conforms to in this superclass protein feature, is the th odds ratio of the subclass feature.

[0053] Furthermore, based on the evaluation value of the molecular phenotypic feature and the comprehensive evaluation value of the feature at the superclass protein level, molecular phenotypic feature annotation and superclass protein feature annotation are performed on each missense mutation in the training dataset to obtain a molecular phenotypic feature annotation result and a superclass protein level feature annotation result respectively, specifically including:

[0054] 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 subclass protein level features that the missense mutation conforms to as the evaluation values of the corresponding subclass protein level features, and label the subclass protein level features that do not conform to as 0 to obtain the annotation result of the subclass protein level features, and sum the annotation results of the subclass protein level features corresponding to obtain the annotation result of the major class protein level features.

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

[0056] 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;

[0057] Based on 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;

[0058] Based on 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 class protein level feature, sum the evaluation values of all subclass protein level features that the protein where the missense mutation to be tested is located conforms to in each major class protein level feature to obtain the comprehensive evaluation value of this major class protein level feature, and use the comprehensive evaluation values of all major class protein level features as the protein level feature judgment result;

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

[0060] 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 the functional effect classification result of the missense mutation to be tested.

[0061] 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;

[0062] 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.

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

[0064] The present invention combines the evaluation values of 27 molecular phenotypic features, the comprehensive evaluation values of each of the 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 solid 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 the field of precision medicine such as drug target screening, clinical diagnosis, and personalized treatment. BRIEF DESCRIPTION OF THE DRAWINGS

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

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

[0067] Figure 3 It is a schematic diagram of the comparison of the ROC curve (Receiver Operating Characteristic curve) test results between the missense mutation functional effect prediction model in an embodiment of the present invention and other models on the pathogenic mutation independent test data set;

[0068] Figure 4 It is a schematic diagram of the comparison of the ROC curve test results between the missense mutation functional effect prediction model in an embodiment of the present invention and other models on the cancer driver mutation independent test data set;

[0069] 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

[0070] 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.

[0071] Embodiment 1:

[0072] 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:

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

[0074] In this embodiment, missense mutation data was 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 the 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).

[0075] Preliminary screening was performed on the collected ClinVar mutations and Humsavar mutations: only ClinVar missense mutations with clinical significance of "pathogenic / likely pathogenic" or "benign / likely benign" were retained, and mutations with conflicting clinical significance were excluded; only Humsavar missense mutations labeled as "likely pathogenic or pathogenic" and "likely benign or benign" were retained.

[0076] 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 were classified as pathogenic mutations; the "likely benign or benign" mutations in ClinVar and the "likely benign or benign" mutations in Humsavar were classified as benign mutations.

[0077] 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 the present 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:

[0078] For mutations with the original identifier being a 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.

[0079] 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.

[0080] For mutations with original identifiers as 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.

[0081] After completing the mutation standardization, observe whether there are identical UniProt canonical mutations that occur 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 the mutation in the dataset, and only retain the repeated UniProt canonical mutation once.

[0082] For the mutations in the dataset, according to subsequent statistical analysis, 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. Therefore, proteins with sequence lengths exceeding 3,000 were removed to avoid potential biases introduced by overly long proteins in subsequent calculations while retaining the vast majority of proteins.

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

[0084] 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 few 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.

[0085] 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.

[0086] 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.

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

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

[0089] 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.

[0090] Therefore, the 1st to 2nd 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.

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

[0092] Therefore, the 3rd to 5th 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.

[0093] 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 to the AlphaFold2 structure and use Python to calculate the shortest spatial distance between the beta-carbon atom (if it is glycine, the alpha-carbon atom is used) of each non-PTM residue in the structure and the PTM site based on the residue atomic coordinates of 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.

[0094] Therefore, the 6th to 8th 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.

[0095] 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 9th residue-level feature is that the mutation occurs at a disordered site.

[0096] 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 the mutation occurred at a short linear motif site.

[0097] Allosteric site: Experimentally confirmed allosteric sites in the Protein Data Bank (PDB) were obtained from the Allosteric Database (ASD). 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 standard sequences.

[0098] 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 the 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 standard 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 standard protein sequences. The allosteric site data from these three sources were integrated to obtain 1,435,333 allosteric site feature data of 14,160 proteins. The 11th residue-level feature was whether the mutation occurred at an allosteric site.

[0099] Protein-protein binding sites: Experimental and predicted protein-protein binding site data were collected or calculated from the following two sources: (1) Experimental 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 homo-dimer structures of 3,166 proteins and 2,445 hetero-dimer 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 for protein interaction interfaces), and then the structures and protein standard sequences were mapped through sequence alignment to obtain 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 the interface residues specific to all interaction partners of each protein were integrated as protein-protein binding sites to obtain 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.

[0100] 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 occurred. 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 protein standard sequence, resulting in 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.

[0101] 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 mutated chains were also considered. A total of 3,714 mutated chains interacting with nucleic acids were obtained, involving 394 proteins. Through the EBI-PISA API, the protein-nucleic acid binding sites in the mutated chains were extracted, and the structural sites were mapped to the sequences, obtaining 15,526 nucleic acid binding sites for 394 proteins. (2) Predicted data from GraphBind (a computational tool name): 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, based on the subsequent results, only the experimental data were retained as the final protein-nucleic acid binding site feature data. Therefore, the 14th residue-level feature is that the mutation occurs at the protein-nucleic acid binding site.

[0102] 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.

[0103] Physical and chemical properties of amino acid residues (6 groups): According to the basic physical and chemical 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.

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

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

[0106] Change in amino acid residue charge before and after mutation (3 groups): Charge values are assigned to 20 standard amino acid residues. Among them, 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 amino acid residue charge 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).

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

[0108] Change in amino acid residue volume before and after mutation (3 groups): The 20 standard amino acid residues are divided into three groups according to 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 amino acid residue volume 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).

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

[0110] 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 (OR) is a statistical metric used to evaluate the likelihood of an event occurring in one group relative to another group. The association difference between pathogenic mutations and benign mutations in the 21 residue-level features and 6 mutation-level features is quantified through the odds ratio. The data calculation expression for OR is:

[0111]

[0112] where, represents the odds ratio of pathogenic missense mutations and benign missense mutations in a molecular phenotype feature, represents the number of pathogenic missense mutations that conform to this molecular phenotype feature, represents the number of pathogenic missense mutations that do not conform to this molecular phenotype feature, represents the number of benign missense mutations that conform to this molecular phenotype feature, represents the number of benign missense mutations that do not conform to this molecular phenotype feature.

[0113] The 27 molecular phenotype 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, three changes in residue volume before and after mutation, and features based on structure, such as core, surface, three secondary structures, PTM_4Å, PTM_8Å, allosteric sites, protein-protein binding sites, protein-small molecule ligand binding sites, protein-nucleic acid binding sites.

[0114] 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 structure features obtained from AlphaFold2 predicted structures or crystal structures, 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 it should be noted whether the standard sequence region covered by the AlphaFold2 structure has PTM_4Å residues to determine whether the protein is included in the analysis.

[0115] 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. , , and These four statistics are used as inputs, and the odds ratio (OR) and p-value are calculated using the fisher_exact() function based on the two-sided Fisher's exact test method.

[0116] To control the influence of multiple hypothesis testing, the Benjamini-Hochberg (BH) multiple testing correction method was used to correct the p-value obtained from Fisher's exact test, resulting in the corrected q-value. The 95% confidence interval was 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, factors such as sample size and effect size. Generally, if the OR value is not equal to 1, it indicates that the 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 the feature can significantly distinguish pathogenic and benign mutations. Among them, OR > 1 indicates that the feature is significantly associated with pathogenic mutations; OR < 1 indicates that the feature is significantly associated with benign mutations.

[0117] 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 the 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Å from it, and PTM_8Å represents the residue at the post - translational modification site and the residues within a spatial distance of 8Å from it.

[0118] 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.

[0119] 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 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.

[0120] 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 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

[0121] indicating a decrease in their binding affinity, and the protein - protein interaction is disrupted to a certain extent. Figure 5 Figure 5 ​Protein stability changes caused by pathogenic missense mutations (Pathogenic) and benign missense mutations (Benign) ( ), and binding affinity changes ( ), Figure 5 In a and b, to measure the or 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, then 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.

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

[0123] 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.

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

[0125] Subcellular localization: The subcellular localization of human proteins 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, a total of 36 subcellular localization information for 11,926 proteins was collected. In this example, the 36 subcellular localization information is referred to as 36 subclass features, and it was respectively determined and analyzed whether the protein where the mutation occurred conforms to these 36 subclass features.

[0126] 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 referred to as 186 subclass features, and it was respectively determined and analyzed whether the protein where the mutation occurred conforms to these 186 subclass features.

[0127] 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 referred to as 32 subclass features, and it was respectively determined and analyzed whether the protein where the mutation occurred conforms to these 32 subclass features.

[0128] Protein specific domains: A total of 4,949 types of specific protein domains of humans were obtained from InterPro. In this example, the 4,949 types of specific domain types are referred to as 4,949 subclass features, and it was respectively determined and analyzed whether the protein where the mutation occurred conforms to these 4,949 subclass features.

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

[0130]

[0131] Where represents the odds ratio of pathogenic missense mutations and benign missense mutations in a subclass of protein level features, represents the number of pathogenic missense mutations that conform to this subclass of protein level features, The number of pathogenic missense mutations that do not conform to the protein-level characteristics of this subclass The number of benign missense mutations that conform to the protein-level characteristics of this subclass The number of benign missense mutations that do not conform to the protein-level characteristics of this subclass.

[0132] 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.

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

[0134]

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

[0136] 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 categories of protein-level characteristics, and the predicted values of the changes in protein stability and binding affinity caused by the mutations are used to form the feature matrix of the model training dataset.

[0137] Finally, a random forest classification model was constructed, with the feature matrix of the model training dataset as input. The optimal hyperparameter combination was determined through a systematic hyperparameter optimization strategy to maximize the model's predictive performance. Specifically, two key hyperparameters were 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, ranging from 50 to 1950, with a step size of 50. This range covers forest configurations from small to large scales to fully evaluate the impact of different numbers of trees on model performance. Max_features controls the maximum number of features considered during the construction of each tree, and its values ​​include strategies of using all features (None), selecting features by the square root of the feature (sqrt), and selecting features by logarithm (log2). In order to fully explore the configuration space of these hyperparameters, a parameter grid containing all possible parameter combinations was constructed, and a five-fold cross-validation of the training set was performed on this basis. This strategy ensures that the impact of each combination on model performance can be systematically evaluated to select the optimal hyperparameter combination. Finally, when the number of decision trees was equal to 1950 and the strategy of logarithmic feature selection was used, the best five-fold cross-validation average Matthews Correlation Coefficient (MCC) was obtained, and the missense mutation functional effect prediction model MutaPheno was obtained.

[0138] Next, we will use an independent test set to evaluate the performance of the missense mutation functional effect prediction model and other existing models.

[0139] In order to compare the prediction effects of different models, the following evaluation indicators are introduced: Recall, also known as True Positive Rate (TPR), which is used to indicate the proportion of positive mutations that the model can correctly identify, that is, the proportion of positive mutations that are predicted to be positive. Its expression is:

[0140] ;

[0141] in, represents the recall rate, represents the true positives, that is, the number of mutations that are actually positive and predicted to be positive, represents the false negatives, i.e., the number of mutations that are actually positive but predicted to be negative.

[0142] The false positive rate (FPR) indicates the proportion of negative mutations that are misclassified as positive mutations. Its expression is:

[0143] ;

[0144] 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.

[0145] Precision, which represents how many of the mutations predicted by the model as positive are truly positive mutations, that is, the proportion of positive mutations among the mutations predicted by the model as positive. Its expression is:

[0146] ;

[0147] Among them, 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.

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

[0149] ;

[0150] Among them, represents the F1 score value, represents recall, represents precision.

[0151] MCC (Matthews correlation coefficient), which is also an index for comprehensively evaluating the performance of the model and is reliable in the case of imbalanced positive and negative data. Its expression is:

[0152] ;

[0153] Among them, 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.

[0154] 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 a model to distinguish between positive and negative samples. AUROC is the area under the ROC curve (Receiver Operating Characteristic curve). 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.

[0155] 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 the two sources, and the operations were the same as those described above, so they will not be elaborated here. The standardized mutation data were merged, and only the missense mutations that were successfully mapped to the structure were retained, and the mutations that were repeated with the training set were removed to ensure that the test set was independent. 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.

[0156] 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.

[0157] To further ensure the fairness of the evaluation and eliminate the bias caused by differences in mutation predictability between 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.

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

[0159]

[0160] It is worth noting that 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. The ROC curve plots the relationship between TPR and FPR.

[0161] 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.

[0162] 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).

[0163] OncoKB is an evidence-based precision oncology knowledge base that provides researchers with detailed somatic mutation annotation information, including the biological effects and cancer-driving nature of mutations, and retains missense mutations annotated as "oncogenic / likely oncogenic" or "likely neutral". The identifier provided by OncoKB mutations is the transcript of the genomic database (Ensembl).

[0164] CGI is a multi-functional platform for evaluating and grading the driver effects of tumor mutations based on existing knowledge and evidence. It integrates data from multiple databases and literature sources. The identifier provided by CGI mutations is the Ensembl transcript.

[0165] The comprehensive cancer missense mutation dataset is integrated from multiple sources and annotated based on the experimental verification 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 identifier provided by the mutations is the gene symbol.

[0166] 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.

[0167] Mutation standardization and protein structure mapping are performed on the mutation data from the three sources in the same way as above, which will not be elaborated here. After merging the standardized mutations, only the missense mutations successfully mapped to the structure are retained, and the mutations 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.

[0168] 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.

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

[0170] 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 model MVP that showed the best performance in predicting driver mutations in the dbNSFP tool; 4 tools specifically for predicting cancer driver mutations (CHASM, CHASMplus, ParsSNP, AI-Driver); and MutaPheno proposed by the present invention.

[0171] 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.

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

[0173]

[0174] MutaPheno achieved the highest area under the receiver operating characteristic curve (AUROC) on this dataset, which was 0.841, and maintained the overall first rank, 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.

[0175] 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.

[0176] Example 2:

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

[0178] 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.

[0179] 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 small category 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 small categories 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 small categories, and obtain the five major categories of protein-level characteristic judgment results.

[0180] 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.

[0181] 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 input, and based on the missense mutation functional effect prediction model MutaPheno described in Example 1, perform functional effect prediction on the missense mutation to be tested, and output the functional effect classification label of the missense mutation to be tested:

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

[0183] The above is only the preferred implementation manner of the present invention. It should be pointed out 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 missense mutation functional effect prediction model, characterized in that Comprising: Integrate missense mutation data from different sources to construct a training dataset, where the missense mutations include pathogenic missense mutations and benign missense mutations; Perform enrichment analysis on the molecular phenotype features and protein-level features of each missense mutation in the training dataset to respectively 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; 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 respectively obtain the molecular phenotype feature annotation result and the major protein-level feature annotation result; 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; Construct a random forest classification model; Use the multi-dimensional feature vector as input to train the random forest classification model 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, wherein 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, wherein After constructing the training dataset, it also includes performing missense mutation standardization and protein structure mapping on the training dataset; 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: 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; 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; When the identifier of the missense mutation is a transcript of the reference sequence database, a transcript of the genomic database, or a protein of the genomic database, convert the transcript of the reference sequence database, the transcript of the genomic database, or the protein of the genomic database where the missense mutation is located into the isoform 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 isoform protein sequence of the unified protein resource database. If they are consistent, retain them; if not, delete the mutation; Perform a pairwise sequence alignment between the standard protein sequence of the UniProt database and the isoform protein sequence of the UniProt database. Map the mutation positions originally based on the isoform protein sequence of the UniProt database to the mutation positions based on the standard protein sequence of the UniProt database, and check whether the reference residues of the missense mutations are consistent with the residues at the corresponding missense mutation positions in the standard protein sequence of the UniProt database. If they are consistent, retain them; if not, delete the missense mutations. The retained missense mutations are used as the standardized missense mutations. The protein structure mapping includes: Map the missense mutations to the AlphaFold2 structure sequence according to the positions and reference residues of the standardized missense mutations in the standard protein sequence of the UniProt database. Among them, if the standard protein sequence of the UniProt database of the protein where the missense mutation is located is the same as the AlphaFold2 structure sequence, map directly; if the standard protein sequence of the UniProt database of the protein where the missense mutation is located is different from the AlphaFold2 structure sequence, perform a pairwise sequence alignment between the standard protein sequence of the UniProt database and the AlphaFold2 structure 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 for each missense mutation data in the training dataset to obtain the molecular phenotype characteristic evaluation value includes: 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. Use 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. Use the logarithm function with base 10 of the odds ratio of the molecular phenotype characteristic as the evaluation value of this molecular phenotype characteristic.

5. The method for constructing a missense mutation functional effect prediction model according to claim 4, characterized in that The molecular phenotype 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 helix type, the secondary structure of the protein where the mutation occurs is a sheet 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 mutation occurs at a residue within a 4 Å spatial distance from the post-translational modification site, the mutation occurs at a residue within an 8 Å spatial distance 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 property of the mutant reference residue belongs to aliphatic, the physicochemical property of the mutant reference residue belongs to aromatic, the physicochemical property of the mutant reference residue belongs to positively charged, the physicochemical property of the mutant reference residue belongs to negatively charged, the physicochemical property of the mutant reference residue belongs to neutral, and the physicochemical property of the mutant reference residue belongs to special; 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, and the absolute value of the volume change of the amino acid residue before and after the mutation is 2.

6. The method for constructing a missense mutation functional effect prediction model according to claim 5, wherein The physicochemical properties of the mutant reference residue are classified based on 20 standard amino acids, specifically including: 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; The method for judging the absolute value of the charge change of the amino acid residue before and after the mutation includes: Setting charge parameters, specifically including: lysine and arginine are +1 positive charges, aspartic acid and glutamic acid are -1 negative charges, and the remaining amino acids are 0 neutral charges; Calculating the absolute value of the charge change of the amino acid residue before and after the mutation according to the charge parameters and the amino acid types of the amino acid residue before and after the mutation; The method for judging the absolute value of the volume change of the amino acid residue before and after the mutation includes: Setting volume codes, specifically including: the volume codes of alanine, glycine, serine, cysteine, proline, threonine, aspartic acid, and asparagine are 0, the volume codes of valine, histidine, glutamic acid, and glutamine are 1, and the volume codes of isoleucine, leucine, methionine, lysine, arginine, phenylalanine, tryptophan, and tyrosine are 2; Calculating the absolute value of the volume change of the amino acid residue before and after the mutation according to the volume codes and the amino acid types of the amino acid residue before and after the mutation.

7. The method for constructing a missense mutation functional effect prediction model according to claim 1, wherein The above-mentioned large-category protein level features include human protein function types, human protein subcellular localization information types, the standard protein range 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 above-mentioned large-category protein level features includes multiple small-category protein level features, and the comprehensive evaluation value of the large-category protein level feature is obtained by summing up the evaluation values of all small-category protein level features that the protein where the missense mutation is located in each large-category protein level feature conforms to.

8. The method for constructing a missense mutation functional effect prediction model according to claim 7, wherein The enrichment analysis of protein level features for each missense mutation in the training dataset to obtain the comprehensive evaluation value of the large-category protein level feature includes: Calculating the number of pathogenic missense mutations and benign missense mutations that conform to the small-category protein level feature and the number of pathogenic missense mutations and benign missense mutations that do not conform to the small-category protein level feature in the training dataset, and using them as inputs to calculate the odds ratio of each small-category protein level feature based on the two-sided Fisher's exact test method; For any missense mutation, sum the logarithms to the base 10 of the odds ratios of all small-category protein level features that the protein where the missense mutation is located conforms to in each large-category protein level feature, and use it 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 function effect prediction model according to claim 4 or 8, characterized in that, The calculation expression of the odds ratio of the missense mutation in the molecular phenotype feature is: ； 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; The calculation expression of the odds ratio of the missense mutation in the small-category protein level feature is: ； Among them, represents the odds ratio of pathogenic missense mutations and benign missense mutations in a small class of protein-level features, represents the number of pathogenic missense mutations that conform to the protein-level features of this small class, represents the number of pathogenic missense mutations that do not conform to the protein-level features of this small class, represents the number of benign missense mutations that conform to the protein-level features of this small class, represents the number of benign missense mutations that do not conform to the protein-level features of this small class; The calculation expression of the comprehensive evaluation value of the large-category protein level feature is: ； Among them, represents any one of the large-category protein level features, represents the comprehensive value of the odds ratio of pathogenic missense mutations and benign missense mutations in this large-category protein level feature, is the number of all small-category features that the protein where the mutation is located conforms to in this large-category protein feature, is the odds ratio of the small-category feature.

10. The method for constructing a missense mutation functional effect prediction model according to claim 1, wherein Based on the molecular phenotype feature evaluation value and the comprehensive evaluation value of the large-category protein level feature, performing molecular phenotype feature annotation and large-category protein feature annotation on each missense mutation in the training dataset respectively to obtain a molecular phenotype feature annotation result and a large-category protein level feature annotation result, specifically including: 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 a molecular phenotype annotation result; label the small-category protein level features that the missense mutation conforms to as the evaluation values of the corresponding small-category protein level features, and label the small-category protein level features that do not conform to as 0 to obtain a small-category protein level feature annotation result, and sum the corresponding small-category protein level feature annotation results to obtain a large-category protein level feature annotation result.

11. A method for predicting the functional effect of missense mutations, characterized in that, Including: 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; Based on 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 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 the molecular phenotype characteristic judgment result is obtained; Based on 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 protein-level characteristic, sum up the evaluation values of all the minor protein-level characteristics that the protein where the missense mutation to be tested is located conforms to in each major protein-level characteristic to obtain the comprehensive evaluation value of this major protein-level characteristic, and use the comprehensive evaluation values of all major protein-level characteristics as the protein-level characteristic judgment result; Calculate the predicted values of the changes in protein stability and binding affinity caused by the missense mutation to be tested; Use the molecular phenotype characteristic judgment result, the protein-level characteristic judgment result, and the predicted values of the changes in protein stability and binding affinity as inputs, and perform functional effect prediction on the missense mutation to be tested using the missense mutation functional effect prediction model obtained based on the missense mutation functional effect prediction model construction method described in any one of claims 1 to 10, and output the functional effect classification result of the missense mutation to be tested.

12. The method for predicting the functional effect of a missense mutation according to claim 11, wherein, 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 driver missense mutation; the missense mutation without functional impact includes benign missense mutation or passenger missense mutation.

Images