Gene prediction method based on co-evolutionary feature analysis

Abstract

The present invention provides a crop genome gene prediction method based on co-evolution feature analysis. The method constructs a multi-species homologous sequence data set, calculates weighted mutual information values, establishes a co-evolution feature network, and combines network topology features and sequence conservation features to use a deep learning model to predict potential coding regions. The method of the present invention integrates sequence alignment quality scoring, adaptive network pruning and multi-level feature analysis, thereby improving the accuracy and reliability of gene prediction in crop genomes. The method can be used to identify promoter regions, exon regions and terminator regions, providing an effective tool for crop genome annotation.

Description

Technical Field

[0001] The invention belongs to the technical field of gene analysis, and in particular relates to a gene prediction method based on co-evolution feature analysis. Background Art

[0002] Accurately identifying and predicting the functional regions of genes, especially coding regions such as promoters, exons and terminators, is of great significance for understanding the regulatory mechanisms and functions of genes. However, traditional gene prediction methods mainly rely on sequence conservation characteristics and fail to fully utilize the co-evolutionary information between sequences, resulting in limited prediction accuracy.

[0003] Current gene prediction methods mostly rely on primary structural information of sequences or simple conservation analysis, ignoring the co-evolutionary relationships between sequence sites. These co-evolutionary relationships often contain important functional information. When performing multiple sequence alignments, traditional methods rarely consider the quality differences of sequence alignments, which may lead to the introduction of noise information and affect prediction accuracy. When using existing methods to predict crop genomes, they often cannot produce good results. For example, for the rice genome, the rice genome has a complex structure and contains a large number of repetitive sequences and transposon elements. There are significant genomic differences between rice subspecies, and it is necessary to consider the impact of population genetic diversity on gene structure. Rice gene expression is tissue-specific and spatiotemporal-specific. These regulatory characteristics are often reflected in subtle changes in sequences, which increases the difficulty of accurate prediction.

[0004] Therefore, there is an urgent need for a gene prediction method that can comprehensively utilize sequence co-evolution characteristics and conservative characteristics to improve the accuracy of gene prediction. Summary of the invention

[0005] The present invention provides a gene prediction method based on co-evolution feature analysis, combines network topology features and sequence conservation features, and uses a deep learning model to achieve the purpose of predicting potential coding regions.

[0006] The scheme of the present invention is as follows:

[0007] A gene prediction method based on co-evolution feature analysis, the method comprising the following steps:

[0008] Step S1, constructing a multi-species homologous sequence dataset, and using a multiple sequence alignment algorithm to align the homologous sequence dataset to obtain an alignment quality score.

[0009] Step S2, calculating the mutual information value between sequence sites, obtaining the weighted mutual information value based on the alignment quality score, constructing a co-evolution feature network according to the significance analysis result of the weighted mutual information value, and performing network pruning.

[0010] Step S3, based on the pruned co-evolution feature network, extract the network topology features, calculate the sequence conservation features, and combine the two types of features to form a comprehensive feature vector.

[0011] Step S4, using a deep learning model to build a classifier, and using the comprehensive feature vector to train classifier parameters.

[0012] Step S5, extracting the feature vector of the sequence to be predicted and inputting it into the classifier, outputting the position information of the potential coding region, wherein the coding region includes the promoter region, the exon region and the terminator region.

[0013] Furthermore, the specific steps of using a multiple sequence alignment algorithm to align the homologous sequence data set to obtain an alignment quality score are:

[0014] Step S11, calculate the distance matrix between sequences ,in Representation sequence and The edit distance of , is the total number of sequences in the homologous sequence dataset.

[0015] Step S12, based on the distance matrix , UPGMA clustering method was used to construct the guide tree .

[0016] Step S13, according to the guide tree The sequences or sequence groups are aligned in sequence according to the branching order.

[0017] Step S14, using the weighted sum formula to calculate the comparison quality score: ,in For the The weight coefficient of a sequence is obtained by normalizing the sequence similarity. For the The score of a sequence ranges from [0 to 1] and is calculated using the BLOSUM62 substitution matrix.

[0018] Further, the calculating the mutual information value between sequence sites and obtaining the weighted mutual information value based on the comparison quality score comprises:

[0019] Statistical sites and The joint frequency distribution of For site and The combined frequency of base pairs above, where ,symbol" ” indicates an empty space.

[0020] Then the marginal frequency distribution is: , ; For site The marginal probability distribution of the bases, For site The marginal probability distribution of the bases, is the sum variable, indicating a traversal of all possible base types.

[0021] Location and The mutual information value between: ; Combined with the comparison quality score, calculate the weighted mutual information value: ,in Site pair The corresponding alignment quality score.

[0022] The normalized mutual information value is obtained by Z-score normalization:

[0023] ,in, is the mean of the weighted mutual information values, is the standard deviation of the weighted mutual information value.

[0024] Furthermore, the constructing of a co-evolution feature network according to the significance analysis results of the weighted mutual information value and performing network pruning includes:

[0025] The normalized mutual information value As a site in a coevolutionary feature network and The significance score of the edge formed by the connection, setting the adaptive threshold ,when Delete the edge when is the benchmark threshold, with a value range of [0.5, 2], and is adaptively adjusted according to the characteristics of the dataset.

[0026] Furthermore, the network topology features include: degree centrality of nodes, weighted betweenness centrality of nodes, weighted closeness centrality of nodes, network clustering coefficient considering mutual information weights, and average shortest path length of the network considering mutual information weights.

[0027] The degree centrality of the node is used to characterize the co-evolutionary importance of the site: ,in, For Node The degree of N is the total number of network nodes.

[0028] The weighted betweenness centrality of the node is used to characterize the mediating role of the site in the co-evolutionary network: ,in, Passing Node Time Node and The number of shortest paths between For Node and The total number of shortest paths between and .

[0029] The weighted proximity centrality of the node is used to characterize the degree of co-evolution of the site with other sites: ,in, For Node To Node The shortest distance.

[0030] The network clustering coefficient considering the mutual information weight is used to characterize the compactness of the local co-evolutionary module:

[0031] ;in, For Node The actual number of connections between neighbors, For Node The average weighted mutual information value of the correlation.

[0032] The average shortest path length of the network considering the mutual information weight is used to characterize the propagation efficiency of the global co-evolution signal:

[0033] ,in, For Node and The shortest distance between.

[0034] Furthermore, the sequence conservation characteristics include: Shannon entropy considering alignment quality, nucleotide diversity considering alignment quality, and evolutionary rate considering alignment quality.

[0035] The Shannon entropy of the alignment quality is considered as: ; The nucleotide diversity considered for alignment quality: , the evolution rate considering the alignment quality: ,in, is the non-synonymous substitution rate, is the synonymous substitution rate.

[0036] Furthermore, the classifier adopts a deep learning model, and the deep learning model includes: an input layer, a hidden layer and an output layer;

[0037] The input layer receives the feature vector ,in Contains network topology features and sequence conservation features.

[0038] The hidden layer includes three fully connected layers. The calculation formula of the layer is: ,in, is the weight matrix, is the bias vector, , is the input feature vector, and the dimensions of each hidden layer are , V is the input feature dimension, and the initialization of the weight matrix is ​​adapted to the magnitude of the input features.

[0039] The output layer uses the Softmax function to calculate the probability of each category: ,in, No. The weight vector of the class; Indicates The dot product of the class weight vector and the output of the last hidden layer, The value range is {1,2,3}, Indicates the type of coding region, including promoter, exon and terminator, .

[0040] Furthermore, the evaluation indicators of the gene prediction results include: sensitivity, specificity, F1 score and Matthew correlation coefficient.

[0041] Sensitivity: ;

[0042] Specificity: ;

[0043] F1 score: ;

[0044] Matthew correlation coefficient: ;in, ; TP represents the number of true positives, TN represents the number of true negatives, FP represents the number of false positives, and FN represents the number of false negatives.

[0045] Compared with the prior art, the present invention has the following beneficial effects:

[0046] The present invention constructs a multi-species homologous sequence data set and introduces a comparison quality scoring mechanism, combines weighted mutual information value calculation and adaptive network pruning strategy, and establishes a prediction model based on co-evolution characteristics, which not only improves the accuracy of gene prediction in crop genomes, but also enhances the adaptability of the prediction method to sequence data of different qualities; by integrating the comprehensive feature vector formed by network topology features and sequence conservation features, combined with a specially designed deep learning classification architecture, the method can simultaneously capture the local and global evolutionary patterns of the sequence, significantly improving the reliability of the prediction results. BRIEF DESCRIPTION OF THE DRAWINGS

[0047] Figure 1 It is a flow chart of the gene prediction method based on co-evolution feature analysis of the present invention;

[0048] Figure 2 It is a schematic diagram of the prediction results of the rice OsMADS1 gene of the present invention. DETAILED DESCRIPTION

[0049] In order to make the purpose, technical solution and advantages of the present invention clearer, the technical solution of the present invention is described clearly and completely below. Obviously, the described embodiments are part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0050] like Figure 1 As shown, the gene prediction method based on co-evolution feature analysis of the present invention comprises the following steps:

[0051] Step S1, constructing a multi-species homologous sequence dataset, and using a multiple sequence alignment algorithm to align the homologous sequence dataset to obtain an alignment quality score.

[0052] The specific steps of using a multiple sequence alignment algorithm to align the homologous sequence data set to obtain an alignment quality score are:

[0053] Step S11, calculate the distance matrix between sequences ,in Representation sequence and The edit distance of , is the total number of sequences in the homologous sequence dataset.

[0054] Step S12, based on the distance matrix , UPGMA clustering method was used to construct the guide tree .

[0055] Step S13, according to the guide tree The sequences or sequence groups are aligned in sequence according to the branching order.

[0056] Step S14, using the weighted sum formula to calculate the comparison quality score: ,in For the The weight coefficient of a sequence is obtained by normalizing the sequence similarity. For the The score of a sequence ranges from [0 to 1] and is calculated using the BLOSUM62 substitution matrix.

[0057] In this embodiment, the genome data of 8 species including rice and its closely related species, including common wild rice (Oryza rufipogon), black rice (Oryza officinalis), and medicinal rice (Oryza medicinal) were selected. From these species, 40 known functionally verified genes and their upstream and downstream 2000bp sequences were selected to construct a homologous sequence data set. Multiple sequence alignment was performed using MUSCLE software, and the parameters were set as: -maxiters 16 -diags -sv -distance1kbit20_3.

[0058] In the analysis of the rice OsMADS1 gene family, a total of 32 homologous sequences from 8 species were first collected. When calculating the edit distance between sequences, the improved Needleman-Wunsch algorithm was used, and the distance matrix obtained showed that the similarity between the two closest sequences was 92%, and the similarity between the two farthest sequences was 65%. The guide tree constructed using the UPGMA method showed that sequences from the same species tended to cluster first. In the final alignment quality assessment, the scores of the conserved regions were generally above 0.85, while the scores of the variable regions were mostly between 0.3-0.6. For the promoter region of the OsMADS1 gene, multiple highly conserved cis-acting elements were found, and their alignment quality scores were all above 0.9.

[0059] Step S2, calculating the mutual information value between sequence sites, obtaining the weighted mutual information value based on the alignment quality score, constructing a co-evolution feature network according to the significance analysis result of the weighted mutual information value, and performing network pruning.

[0060] The calculating the mutual information value between sequence sites and obtaining the weighted mutual information value based on the comparison quality score comprises:

[0061] Statistical sites and The joint frequency distribution of For site and The combined frequency of base pairs above, where ,symbol" ” indicates an empty space.

[0062] Then the marginal frequency distribution is: , ; For site The marginal probability distribution of the bases, For site The marginal probability distribution of the bases, is the sum variable, indicating a traversal of all possible base types.

[0063] Location and The mutual information value between: ; Combined with the comparison quality score, calculate the weighted mutual information value: ,in Site pair The corresponding alignment quality score.

[0064] The normalized mutual information value is obtained by Z-score normalization:

[0065] ,in, is the mean of the weighted mutual information values, is the standard deviation of the weighted mutual information value.

[0066] The method of constructing a co-evolution feature network according to the significance analysis result of the weighted mutual information value and performing network pruning includes:

[0067] The normalized mutual information value As a site in a coevolutionary feature network and The significance score of the edge formed by the connection, setting the adaptive threshold ,when Delete the edge when is the benchmark threshold, with a value range of [0.5, 2], and is adaptively adjusted according to the characteristics of the dataset.

[0068] In the process of constructing the rice gene network, the optimal benchmark threshold β=1.5 was determined by analyzing the known 40 gene data sets. When using this threshold for network pruning, about 28% of the edges in the OsMADS1 gene were retained, and these edges were mainly distributed in the known functional regions. In particular, within the ±50bp range near the transcription start site, the density of retained edges was 2.3 times that of other regions. By comparing with ChIP-seq data, it was found that the nodes connected by the retained edges had significant overlap with the transcription factor binding sites, verifying the rationality of the pruning strategy.

[0069] In the analysis of rice gene sequences, it was found through calculation that the mutual information value between site pairs in known exon regions is significantly higher than that in other regions. For example, in the OsMADS1 gene, the average mutual information value of the exon region is 0.72, while the average value of the intron region is only 0.31. When pruning the network, the benchmark threshold β=1.2 is set, which is the optimal parameter determined by cross-validation of known gene sets. About 35% of the edges are retained in the pruned network, which are mainly concentrated in the regulatory elements and coding regions of functional verification.

[0070] Taking the first exon region of rice OsMADS1 gene as an example, the mutual information values ​​of all site pairs in this region were analyzed. The statistical results showed that about 15% of the site pairs in this region showed significant co-evolutionary relationships (normalized mutual information values ​​greater than 2). Especially in the functional domain encoding region, such as the MADS-box domain, the mutual information values ​​between sites were significantly higher than those in other regions. For example, within the MADS-box domain, the average weighted mutual information value was 0.82, while in the non-domain region, the value was only 0.35. This difference also has similar manifestations in other known functional regions, providing an important basis for predicting new functional elements.

[0071] Step S3, based on the pruned co-evolution feature network, extract the network topology features, calculate the sequence conservation features, and combine the two types of features to form a comprehensive feature vector.

[0072] The network topology features include: degree centrality of nodes, weighted betweenness centrality of nodes, weighted closeness centrality of nodes, network clustering coefficient considering mutual information weight, and average shortest path length of the network considering mutual information weight.

[0073] The degree centrality of the node is used to characterize the co-evolutionary importance of the site: ,in, For Node The degree of N is the total number of network nodes.

[0074] The weighted betweenness centrality of the node is used to characterize the mediating role of the site in the co-evolutionary network: ,in, Passing Node Time Node and The number of shortest paths between For Node and The total number of shortest paths between and .

[0075] The weighted proximity centrality of the node is used to characterize the degree of co-evolution of the site with other sites: ,in, For Node To Node The shortest distance.

[0076] The network clustering coefficient considering the mutual information weight is used to characterize the compactness of the local co-evolutionary module:

[0077] ;in, For Node The actual number of connections between neighbors, For Node The average weighted mutual information value of the correlation.

[0078] The average shortest path length of the network considering the mutual information weight is used to characterize the propagation efficiency of the global co-evolution signal:

[0079] ,in, For Node and The shortest distance between.

[0080] The network topology analysis of the rice OsMADS1 gene showed that the network characteristics of the nodes in the known key regulatory regions have significant characteristics. For example, in the TATA-box region upstream of the transcription start site, the average degree centrality is 0.68, which is significantly higher than the 0.23 of the random region; the weighted betweenness centrality of this region also reaches 0.55, indicating that it plays an important mediating role in the network.

[0081] The sequence conservation characteristics include: Shannon entropy considering alignment quality, nucleotide diversity considering alignment quality, and evolution rate considering alignment quality.

[0082] The Shannon entropy of the alignment quality is considered as: ; The nucleotide diversity considered for alignment quality: , the evolution rate considering the alignment quality: ,in, is the non-synonymous substitution rate, is the synonymous substitution rate.

[0083] When analyzing the conservative features in the rice genome, 40 functionally verified genes were selected as references. Among these genes, the average Shannon entropy of the promoter region was 0.75, that of the exon region was 0.85, and that of the non-functional region was only 0.32. Nucleotide diversity analysis considering the alignment quality showed that the nucleotide diversity at key regulatory element locations, such as transcription factor binding sites, was significantly lower than that in the surrounding regions, with an average of 0.15. Evolutionary rate analysis showed that in regions with strong functional constraints, such as near essential splicing sites, the dN / dS ratio was generally less than 0.1, indicating that these positions were subject to strong purifying selection.

[0084] In the analysis of the rice gene OsMADS1, the network topology characteristics of its coding region showed that the average degree centrality was 0.45, the average weighted betweenness centrality was 0.38, and the average weighted proximity centrality was 0.62. These values ​​were significantly higher than those in the non-coding region (0.21, 0.15, and 0.33, respectively). Sequence conservation analysis showed that the Shannon entropy of the gene averaged 0.82 in the promoter region and 0.91 in the exon region, which was highly consistent with the functional region verified by the experiment.

[0085] Step S4, constructing a classifier using a deep learning model, and using the comprehensive feature vector to train classifier parameters. The classifier uses a deep learning model, and the deep learning model includes: an input layer, a hidden layer, and an output layer;

[0086] The input layer receives the feature vector ,in Contains network topology features and sequence conservation features.

[0087] The hidden layer includes three fully connected layers. The calculation formula of the layer is: ,in, is the weight matrix, is the bias vector, , is the input feature vector, and the dimensions of each hidden layer are , V is the input feature dimension, and the initialization of the weight matrix is ​​adapted to the magnitude of the input features.

[0088] The output layer uses the Softmax function to calculate the probability of each category: ,in, No. The weight vector of the class; Indicates The dot product of the class weight vector and the output of the last hidden layer, The value range is {1,2,3}, Indicates the type of coding region, including promoter, exon and terminator, .

[0089] In the practice of rice gene prediction, the input feature dimension V is set to 128, which is the optimal parameter determined by grid search. The dimensions of the three hidden layers of the model are 256, 128, and 5, respectively. When initializing the weights, considering the distribution characteristics of the input features, the He initialization method is used so that the variance of the initial weights is 2 / n, where n is the input dimension. In the actual training process, the Adam optimizer is used, the initial learning rate is set to 0.001, and the learning rate decay strategy is adopted, decaying to 0.8 times the original every 20 epochs. To prevent overfitting, a dropout layer is added after each fully connected layer, and the dropout rate is set to 0.3.

[0090] During the training process, feature data from 40 known rice genes were used as the training set, of which 70% was used for training and 30% for validation. The model training used a batch size of 32, a learning rate of 0.001, and 100 training rounds. The best performance indicators achieved on the validation set were: accuracy 92.3%, recall 89.7%, and F1 value 90.9%. In particular, the prediction accuracy for the exon region reached 94.5%, which shows that the model has a good recognition ability for the gene coding region.

[0091] Step S5, extracting the feature vector of the sequence to be predicted and inputting it into the classifier, outputting the position information of the potential coding region, wherein the coding region includes the promoter region, the exon region and the terminator region.

[0092] This method was applied to predict a 100kb unknown functional region on rice chromosome 4. The prediction results successfully identified five potential gene regions, three of which were consistent with existing EST data, and the other two were verified by RT-PCR experiments to have transcriptional activity. Functional analysis of one of the newly predicted genes showed that the gene may be involved in the regulation of rice drought resistance, and the predicted promoter region showed obvious transcriptional activity enhancement under drought stress.

[0093] The evaluation indicators of the gene prediction results include: sensitivity, specificity, F1 score and Matthew correlation coefficient.

[0094] Sensitivity: ;

[0095] Specificity: ;

[0096] F1 score: ;

[0097] Matthew correlation coefficient: ;in, ; TP represents the number of true positives, TN represents the number of true negatives, FP represents the number of false positives, and FN represents the number of false negatives.

[0098] Taking rice genome data as an example, 100 known genes on chromosome 4 were selected as the test data set for model evaluation. The evaluation results show:

[0099] For promoter region prediction: sensitivity (Se) = 0.89, indicating that 89% of all actual promoter regions were correctly predicted; specificity (Sp) = 0.92, indicating that 92% of all non-promoter regions were correctly identified; F1 score = 0.90, reflecting the overall balance of the prediction; Matthew correlation coefficient (MCC) = 0.86, showing that the prediction results have a strong correlation.

[0100] For exon region prediction: sensitivity (Se) = 0.93, indicating that the model can identify real exons well, specificity (Sp) = 0.95, indicating that the model rarely misclassifies other regions as exons, F1 score = 0.94, showing a high comprehensive prediction ability, and Matthew correlation coefficient (MCC) = 0.91, proving that the prediction results are highly credible.

[0101] For the prediction of terminator region: sensitivity (Se) = 0.87, indicating good recognition ability for terminators, specificity (Sp) = 0.91, indicating that there are few false positive predictions, F1 score = 0.89, showing that the overall prediction effect is satisfactory, and Matthew correlation coefficient (MCC) = 0.85, indicating that the prediction results have good reliability.

[0102] Compared with other mainstream prediction methods, this method has shown obvious advantages in all indicators: compared with traditional methods based on sequence features, the sensitivity of this method has increased by an average of 12%, and compared with existing deep learning methods, the specificity has increased by 8%. In terms of predicting new genes, the F1 score of this method is 0.05 higher than that of existing methods on average. In particular, this method shows obvious advantages when predicting genes with atypical sequence features. For example, when predicting genes with non-traditional splicing sites, the accuracy of this method reaches 85%, while the traditional method is only 70%. This improvement is mainly due to the fact that this method comprehensively considers the co-evolutionary characteristics and conservative characteristics of the sequence.

[0103] like Figure 2As shown, the prediction result of the rice OsMADS1 gene using the method of the present invention is shown, the horizontal axis represents the sequence position (unit: bp), the broken line represents the predicted probability score, and the prediction result shows that: 0-500bp is the promoter region, 500-1500bp and 1500-2000bp are two exon regions, and 2000-2500bp is the terminator region.

[0104] In practical applications, this evaluation system can not only be used to evaluate model performance, but also guide model optimization. For example, by analyzing cases of prediction errors, it was found that the prediction accuracy would slightly decrease in areas with particularly high GC content (>70%). Based on this finding, the model performance can be improved by increasing the proportion of high GC content sequences in the training data. This feedback optimization mechanism enables the model to continuously improve prediction accuracy in practical applications.

[0105] It should be noted that those skilled in the art should understand that various changes and equivalent substitutions may be made to the present invention without departing from the scope of the present invention. In addition, various modifications may be made to the present invention for specific situations or materials without departing from the scope of the present invention. Therefore, the present invention is not limited to the specific embodiments disclosed, but should include all embodiments falling within the scope of the claims of the present invention.

Claims

1. A gene prediction method based on co-evolutionary feature analysis, characterized in that: The method comprises the following steps: Step S1, constructing a multi-species homologous sequence dataset, and using a multiple sequence alignment algorithm to align the homologous sequence dataset to obtain an alignment quality score; Step S2, calculating the mutual information value between sequence sites, obtaining the weighted mutual information value based on the alignment quality score, constructing a co-evolution feature network according to the significance analysis result of the weighted mutual information value, and performing network pruning; Step S3, based on the pruned co-evolution feature network, extract network topology features, calculate sequence conservation features, and combine the two types of features to form a comprehensive feature vector; Step S4, using a deep learning model to build a classifier, and using the comprehensive feature vector to train classifier parameters; Step S5, extracting the feature vector of the sequence to be predicted and inputting it into the classifier, outputting the position information of the potential coding region, wherein the coding region includes the promoter region, the exon region and the terminator region.

2. The gene prediction method according to claim 1, characterized in that The specific steps of using a multiple sequence alignment algorithm to align the homologous sequence data set to obtain an alignment quality score are: Step S11, calculate the distance matrix between sequences ,in Representation sequence and The edit distance of , is the total number of sequences in the homologous sequence dataset; Step S12, based on the distance matrix , UPGMA clustering method was used to construct the guide tree ; Step S13, according to the guide tree The sequence or sequence group is aligned in sequence according to the branching order; Step S14, using the weighted sum formula to calculate the comparison quality score: ,in For the The weight coefficient of a sequence is obtained by normalizing the sequence similarity. For the The score of a sequence ranges from [0 to 1] and is calculated using the BLOSUM62 substitution matrix.

3. The gene prediction method according to claim 2, characterized in that: The calculating the mutual information value between sequence sites and obtaining the weighted mutual information value based on the comparison quality score comprises: Statistical sites and The joint frequency distribution of For site and The combined frequency of base pairs above, where ,symbol" ” indicates an empty space; Then the marginal frequency distribution is: , ; For site The marginal probability distribution of the bases, For site The marginal probability distribution of the bases, is the sum variable, indicating a traversal of all possible base types; Location and The mutual information value between: ; Combined with the comparison quality score, calculate the weighted mutual information value: ,in Site pair The corresponding alignment quality score; The normalized mutual information value is obtained by Z-score normalization: ,in, is the mean of the weighted mutual information values, is the standard deviation of the weighted mutual information value.

4. The gene prediction method according to claim 3, characterized in that: The method of constructing a co-evolution feature network according to the significance analysis result of the weighted mutual information value and performing network pruning includes: The normalized mutual information value As a site in a coevolutionary feature network and The significance score of the edge formed by the connection, setting the adaptive threshold ,when Delete the edge when is the benchmark threshold, with a value range of [0.5, 2], and is adaptively adjusted according to the characteristics of the dataset.

5. The gene prediction method according to claim 4, characterized in that: The network topology features include: degree centrality of nodes, weighted betweenness centrality of nodes, weighted closeness centrality of nodes, network clustering coefficient considering mutual information weight, and average shortest path length of the network considering mutual information weight; The degree centrality of the node is used to characterize the co-evolutionary importance of the site: ,in, For Node The degree of, N is the total number of network nodes; The weighted betweenness centrality of the node is used to characterize the mediating role of the site in the co-evolutionary network: ,in, Passing Node Time Node and The number of shortest paths between For Node and The total number of shortest paths between and ; The weighted proximity centrality of the node is used to characterize the degree of co-evolution of the site with other sites: ,in, For Node To Node The shortest distance; The network clustering coefficient considering the mutual information weight is used to characterize the compactness of the local co-evolution module: ;in, For Node The actual number of connections between neighbors, For Node The average weighted mutual information value of the correlation; The average shortest path length of the network considering the mutual information weight is used to characterize the propagation efficiency of the global co-evolution signal: ,in, For Node and The shortest distance between.

6. The gene prediction method according to claim 5, characterized in that: The sequence conservation characteristics include: Shannon entropy considering alignment quality, nucleotide diversity considering alignment quality, and evolutionary rate considering alignment quality; The Shannon entropy of the alignment quality is considered as: ; The nucleotide diversity considered for alignment quality: , the evolution rate considering the alignment quality: ,in, is the non-synonymous substitution rate, is the synonymous substitution rate.

7. The gene prediction method according to claim 6, characterized in that: The classifier adopts a deep learning model, which includes: an input layer, a hidden layer and an output layer; The input layer receives the feature vector ,in Contains network topology features and sequence conservation features; The hidden layer contains three fully connected layers. The calculation formula of the layer is: ,in, is the weight matrix, is the bias vector, , is the input feature vector, and the dimensions of each hidden layer are , V is the input feature dimension, and the initialization of the weight matrix is ​​adapted to the magnitude of the input features; The output layer uses the Softmax function to calculate the probability of each category: ,in, No. The weight vector of the class; Indicates The dot product of the class weight vector and the output of the last hidden layer, The value range is {1,2,3}, Indicates the type of coding region, including promoter, exon and terminator, .

8. The gene prediction method according to claim 7, characterized in that: The evaluation indicators of the gene prediction results include: sensitivity, specificity, F1 score and Matthew correlation coefficient; Sensitivity: ; Specificity: ; F1 score: ; Matthew correlation coefficient: ;in, ; TP represents the number of true positives, TN represents the number of true negatives, FP represents the number of false positives, and FN represents the number of false negatives.

Images