Enzyme kinetic parameter prediction and virtual screening method based on multi-source feature fusion

Abstract

The invention discloses an enzyme kinetic parameter prediction and virtual screening method based on multi-source feature fusion, which belongs to the field of biological enzyme engineering and biological catalysis, and constructs a multi-modal fusion feature vector by combining a protein sequence and structural information of a substrate molecule and respectively utilizing a ProtT5 protein pre-training language model and a MolT5 model. According to the method, the extreme random tree model is used for performing high-precision prediction on the catalytic kinetic parameters of the enzyme-substrate pair. On the basis, systematic residue mutation is carried out on a candidate sequence similar to a target enzyme sequence at a key site, a mutant database is constructed, batch prediction and sorting are carried out on the mutant database, and finally virtual screening of the target enzyme is achieved. The method has remarkable advantages in the aspects of feature extraction, modeling efficiency and screening precision, and is suitable for large-scale enzyme engineering application.

Description

Technical Field

[0001] This invention belongs to the field of bioenzyme engineering and biocatalysis, specifically relating to a method for constructing an enzyme kinetic parameter prediction model based on multi-source feature fusion, for the efficient screening of target enzymes. Background Technology

[0002] In biocatalysis and enzyme engineering research, improving enzyme catalytic efficiency, optimizing substrate affinity, and screening for highly efficient enzymes capable of catalyzing reactions on natural or non-natural substrates are key aspects of directed evolution and industrial applications. Among these, the catalytic constant Kcat and the Michaelis constant Km are the main kinetic parameters for measuring enzyme catalytic activity and substrate binding affinity, and they are instructive for evaluating enzyme performance. Traditional high-throughput experimental screening methods are costly and time-consuming, making them unsuitable for the needs of rapid iterative optimization.

[0003] In recent years, machine learning and deep learning models based on sequence or structural features have made some progress in predicting enzyme kinetic parameters. However, relying solely on a single data source or model architecture often fails to fully integrate protein sequence, physicochemical properties, and substrate information, thus limiting prediction accuracy and generalization ability. Therefore, there is an urgent need for a prediction model that can efficiently integrate multi-dimensional features and stably output high-precision kinetic parameters, and for building a rapid and scalable enzyme sequence screening method based on this model.

[0004] Currently, the field of enzyme virtual screening mainly relies on three technical paths: traditional computation, AI-assisted methods, and traditional experiments. Each method has significant technical shortcomings: In traditional computational methods, structure-based molecular docking and molecular dynamics simulations are highly dependent on the three-dimensional crystal structure of the enzyme and cannot be carried out without structural information. Furthermore, the scoring function is difficult to accurately capture weak interactions such as hydrogen bonds and hydrophobic interactions between the enzyme and substrate, and the dynamic simulation of the flexible regions of proteins will significantly increase the computation time. In ligand screening methods, quantitative structure-activity relationship (QSAR) models require hundreds to thousands of known active ligand data. For rare reaction types or novel non-natural substrates, the lack of data leads to a sharp drop in model prediction accuracy. Pharmacophore model construction is easily affected by the subjective judgment of researchers on key action sites, and structural similarity matching often ignores the actual binding complementarity between the enzyme and ligand. Although near-attack conformational methods simplify the computational complexity of transition states, they lack NAC-based system design strategies and are difficult to accurately guide the structural remodeling of natural enzymes and the optimization of mutant activities.

[0005] While AI-assisted screening methods have improved efficiency through architectures such as Transformer and Graph Neural Networks (GNNs), they still face multiple bottlenecks. Model performance is limited by the quantity and completeness of high-quality enzyme-substrate reaction data. Existing databases contain inconsistencies in the measurement temperature and pH of some Kcat and Km data, resulting in high noise in the training data. Furthermore, data on novel enzyme families or cross-species enzyme sequences are scarce, significantly reducing the model's cross-scenario generalization ability. The "black box" nature of deep learning models leads to a lack of interpretability in screening results. Even if key amino acid residues are located through attention mechanisms, it is difficult to clarify the specific structure-activity relationship between residues and catalytic activity. Training on large-scale datasets and high-throughput screening still rely on high-performance computing clusters, and the hardware threshold limits their widespread application. Traditional experimental screening relies on high-throughput wet experiments, involving complex processes such as gene mutation library construction, host transformation, single-clone selection, and enzyme activity assays. The screening cycle for a single target can take 6-12 months, costing tens to hundreds of thousands of yuan, and can only cover tens of thousands of candidate sequences, far from exhausting the vast enzyme sequence space. Furthermore, issues such as heterologous protein folding abnormalities and buffer component interference can cause some "highly active" candidate enzymes to fail during experimental validation. In addition, all three methods share common drawbacks: insufficient consistency between virtual screening results and experimental validation, with false positive rates generally exceeding 30%; most methods focus only on basic kinetic parameters such as Kcat and Km, neglecting key industrial application indicators such as enzyme thermostability, protease resistance, and substrate specificity, resulting in screened enzymes that are difficult to adapt to actual production scenarios; and all methods struggle to integrate complex environmental factors such as temperature, pH, and solvent effects, further exacerbating the disconnect between screening results and industrial needs. Summary of the Invention

[0006] The purpose of this invention is to overcome the shortcomings of the existing technology and provide a method for predicting enzyme kinetic parameters and virtual screening based on multi-source feature fusion. This method addresses the problems of long screening cycles, high costs, and narrow coverage in traditional high-throughput experiments in the fields of bioenzyme engineering and biocatalysis. By fusing enzyme sequence semantic features, substrate topological fingerprints, and substrate semantic features to construct a composite feature vector, and combining it with an optimized extreme random tree regression model to improve the prediction accuracy of Kcat / Km, and achieving large-scale and efficient virtual screening through pre-screening with mutant instability index, this invention ultimately provides an enzyme screening solution that is accurately adapted to industrial scenarios, promoting the transformation of biocatalysis technology from laboratory research to industrial applications.

[0007] This invention is achieved using the following technical solution: Step 1. Obtain the enzyme sequence and its corresponding kinetic parameters Kcat and Km from public databases, perform data cleaning, redundancy removal and quality control, divide the training set and test set in an 8:2 ratio, and construct subsets of different intervals based on BLAST similarity; Step 2. Use the ProtT5 protein pre-trained model to extract the semantic embedding vector of the enzyme sequence, use the RDKit tool to generate the Morgan fingerprint of the substrate molecule, combine the MolT5 molecular pre-trained model to extract the semantic embedding vector of the substrate, and concatenate the obtained three types of features to form a composite feature vector. Step 3. Train prediction models for Kcat and Km using extreme random tree regression models respectively, and improve model accuracy and generalization performance through cross-validation and hyperparameter tuning; Step 4. Based on the target enzyme sequence and functional site information, a combination of full-coverage random mutagenesis and directed mutagenesis is used to generate mutant combinations, calculate their instability index and remove unstable sequences ≥40. Step 5. Input the candidate mutants into the trained prediction model, predict the Kcat and Km values ​​in batches, calculate the Kcat / Km ratio, and sort them according to the ratio. By default, the top 10 candidate mutant sequences with the highest ratios are output.

[0008] The technical effects of this invention are as follows: This invention employs a multi-source feature fusion strategy, utilizing ProtT5 to extract semantic embedding vectors of enzyme sequences, RDKit to generate substrate Morgan fingerprints, and MolT5 to extract substrate semantic vectors, concatenating them to form a 1×3072-dimensional composite feature. Combined with an extreme random tree regression model optimized through 10-fold cross-validation and hyperparameter tuning, this significantly improves the prediction accuracy of enzyme kinetic parameters Kcat and Km, as well as the generalization ability across different BLAST similarity intervals. In the mutant screening stage, a mutant library is constructed by combining full-coverage random mutation and directed mutation, and stable candidates are screened based on the instability index. Finally, the optimal mutant is output by ranking the Kcat / Km ratio. Compared with traditional high-throughput methods, this not only significantly reduces enzyme screening costs but also enables efficient virtual screening of large-scale enzyme sequences. Attached Figure Description

[0009] Figure 1 This is a schematic diagram illustrating the framework for relating dataset partitioning to prediction results.

[0010] Figure 2 This is a technical roadmap for multi-source feature extraction.

[0011] Figure 3 The bar chart shows the performance evaluation metrics for the Km prediction model.

[0012] Figure 4 This is a bar chart showing the performance evaluation metrics of the Kcat prediction model.

[0013] Figure 5 This is an overview diagram of the overall technical process of enzyme virtual screening in this invention. Detailed Implementation

[0014] The present invention will be further described below with reference to the accompanying drawings and specific embodiments.

[0015] To overcome the limitations of existing technologies, such as limited prediction accuracy and long screening cycles, this invention proposes a method for enzyme kinetic parameter prediction and virtual screening based on multi-source feature fusion. This method combines protein sequence features with substrate molecule features to establish a high-precision kinetic parameter prediction model, and integrates an automated mutant database to achieve rapid optimization of the target enzyme.

[0016] Reference Figure 1 As shown, the dataset processing logic for enzyme kinetic parameter prediction in this invention presents the independent processing paths of the kcat and Km datasets on the left, the data and subsequent prediction stages are connected through the "feature extraction module" in the middle, and the prediction results corresponding to two data partitioning methods, "based on biased dataset random partitioning" and "based on sequence similarity partitioning (layered according to BLAST similarity intervals)," which intuitively demonstrate the impact of different data partitioning strategies on model evaluation and provide visual support for the rational verification of the dataset for subsequent model training.

[0017] Reference Figure 2 As shown, the feature extraction process for enzyme sequences and substrate molecules is described in detail: For enzyme sequences, an n×1024-dimensional embedding matrix is ​​generated using the ProtT5 protein pre-training model, and then averaged and pooled to obtain a 1×1024-dimensional fixed vector; For substrate molecules, a 1×1024-dimensional Morgan structural fingerprint is generated using the RDKit tool, and a 1×1024-dimensional semantic embedding vector is extracted using the MolT5 molecular pre-training model.

[0018] Reference Figure 3 As shown in the bar chart, after 10-fold cross-validation and hyperparameter tuning, the four key performance indicators of the Km prediction model are presented: coefficient of determination R² = 0.77, root mean square error (RMSE) = 0.70, mean absolute error (MAE) = 0.61, and Pearson correlation coefficient (PCC) = 0.59. These quantitative indicators intuitively reflect the accuracy and reliability of the model's prediction of the Km (Miele constant) parameter, providing data support for the model's application in practical screening.

[0019] Reference Figure 4 As shown in the bar chart, the four performance indicators of the optimized Kcat (catalytic constant) prediction model are presented: R² = 0.80, RMSE = 0.74, MAE = 0.63, and PCC = 0.61. (This is achieved by comparing with...) Figure 3The comparison of the Km model indicators reflects the performance differences and overall optimization effect of the dual prediction models, and verifies the applicability of the extreme random tree regression model in the prediction of key parameters of enzyme kinetics, especially highlighting the higher accuracy of Kcat parameter prediction.

[0020] Reference Figure 5 As shown, the method is divided into four core stages: the first stage, "Data Preparation and Quality Control," covers the collection, cleaning, and partitioning of enzyme and substrate data (8:2 training / test set + BLAST stratification); the second stage, "Feature Extraction and Fusion," clarifies the sources of enzyme sequence features (ProtT5 extraction) and substrate features (topological + semantic extraction) and the fusion method of "splicing together into a 3072-dimensional vector"; the third stage, "Model Construction and Mutant Analysis," demonstrates the parallel training of Kcat / Km dual models and the processing flow of mutant "retrieval-mutation-stability prediction-verification"; and the fourth stage, "Batch Prediction and Screening," presents the entire chain of batch feature generation, model call prediction, and output of Top200 candidates sorted by Kcat / Km.

[0021] Data Acquisition First, datasets containing enzyme IDs (EC), UniProtIDs, enzyme sequences, Km, Kcat, substrate names, and substrate SMILES codes were downloaded from the BRENDA and SABIO-RK enzyme databases, respectively. Incomplete entries and entries with negative kinetic parameter values ​​were removed. These were then combined to form the original dataset, which was used for splitting the training and test sets. The Kcat and Km datasets were randomly shuffled and then divided into training and test sets in an 8:2 ratio. A specified number of sub-test sets with similarity to the training set within different ranges (0%–40%, 40%–60%, 60%–80%, and 80%–100%) were randomly selected from the test set using the BLAST algorithm to comprehensively evaluate the model performance. The corresponding prediction models were trained using the Kcat and Km datasets, respectively. Each dataset should contain the enzyme sequence, specified enzyme kinetic parameters, and substrate SMILES codes.

[0022] Feature extraction Protein embedding vectors were generated using the ProtT5 model developed by the Google DeepMind team. This model, based on a transformer architecture and trained with a self-supervised language model, captures the semantic features of protein sequences and outputs an n×1024-dimensional embedding representation. A fixed-dimensional (1×1024) protein feature vector was obtained through average pooling. Substrate structural features were first generated using the open-source molecular chemistry toolkit RDKit, producing a Morgan molecular fingerprint, essentially a structure-dependent one-hot encoding (1×1024 dimensions). This is commonly used to describe the local topological structure information of molecules. Furthermore, the MolT5 model was used to embed the substrate molecule's SMILES expression. MolT5 is a deep generative model based on the T5 language model architecture proposed by Google and pre-trained on the large-scale molecular SMILES dataset. It effectively extracts semantic features at the molecular level and outputs a 1×1024-dimensional continuous vector. The three vectors were concatenated column-wise to obtain a unified composite feature vector (1×3072 dimensions).

[0023] Model training An extreme random tree regression model (based on the Extra Trees regressor in Scikit-learn) is employed. This method is an enhanced ensemble learning model based on random forests, exhibiting higher bias-variance balance. During training, 10-fold cross-validation is combined with an automatic hyperparameter tuning strategy. The improvement in the coefficient of determination (R²) is used as the evaluation metric to optimize both the Kcat and Km prediction models.

[0024] Quick Filtering Process First, a target enzyme is selected. Then, protein sequences that meet the similarity threshold in the Swiss-Prot database are screened using the BLAST algorithm based on sequence similarity and ranked according to overall similarity. Swiss-Prot is a high-quality, manually annotated protein database that provides detailed functional annotations, domain information, and protein sequences.

[0025] Users can customize the total number of enzymes to be screened, and then use the screening system to obtain functional site information for all selected enzymes, including the target enzyme. Mutation combinations are then performed on all enzyme sequences according to the program's default empirical mutation method or custom mutation settings. The instability index is based on the statistical results of the influence of different dipeptides on protein stability, and in this patent, it is implemented using the open-source bioinformatics tool BioPython library. After calculating the instability index of each mutant sequence, mutant sequences with an instability index greater than or equal to 40 are removed, and kinetic parameters are predicted for the obtained enzyme sequence combinations. Then, the Kcat / Km value is calculated for each enzyme sequence, and they are sorted according to the value.

[0026] This invention mainly includes the following steps: 1. Data preprocessing: Obtain enzyme sequences and corresponding Kcat and Km data from public databases, clean and remove redundancy, and divide the data into training and test sets according to the ratio. At the same time, use the BLAST program to construct subsets of different sequence similarity intervals based on sequence similarity.

[0027] 2. Feature Fusion: The enzyme sequence semantic embedding is extracted using the ProtT5 model, the Morgan fingerprint of the substrate is generated using the RDKit tool, and the global semantic features of the substrate are extracted using the MolT5 model. The three types of features are then concatenated into a composite feature vector.

[0028] 3. Model building and optimization: Using extreme random trees as the regression model, prediction models for Kcat and Km are built respectively, and the accuracy and robustness of the models are improved through cross-validation and hyperparameter optimization.

[0029] 4. Construction of mutant library: A strategy combining random and directed mutations is used to introduce mutant combinations at the functional sites of the target enzyme sequence, and the instability index is calculated to eliminate unstable sequences.

[0030] 5. Virtual screening: Input candidate mutant sequences into the prediction model, predict Kcat and Km parameters in batches, calculate and sort the Kcat / Km ratio, and finally select the sequence with the highest ratio as the candidate result.

[0031] Specifically, the following steps are included: First, construct the dataset: Collect enzyme sequences, corresponding kinetic parameters Kcat and Km, and substrate SMILES codes from the BRANDA and SABIO-RK databases. Clean the data according to the following rules: Remove data with Kcat values ​​less than 10. - ³s - ¹or greater than 10 5 s - ¹, Km value less than 10 - ²μM or greater than 10 6 Abnormal entries with a length of μM were removed; entries with invalid substrate SMILES encoding, enzyme sequences shorter than 50 or longer than 1000 amino acid residues were removed; entries lacking key information were also removed. After cleaning, the Kcat and Km datasets were randomly shuffled and divided into training and testing sets in an 8:2 ratio; a specified number of samples were randomly selected from the testing set using the BLAST algorithm to construct four subsets with sequence similarity of 0%~40%, 40%~60%, 60%~80%, and 80%~100%, ensuring that the data in each subset were evenly distributed in the testing set.

[0032] During the model training phase, an extreme random tree regression model was employed, combined with 10-fold cross-validation and gridded hyperparameter search for hyperparameter tuning. The coefficient of determination (R²) was used as the optimization evaluation metric to train the Kcat and Km prediction models respectively. The feature inputs of the two models were completely identical: the enzyme sequence semantic embedding vector generated by the ProtT5 protein pre-training model, the substrate Morgan fingerprint generated by the RDKit tool, and the substrate semantic embedding vector extracted by the MolT5 molecular pre-training model. These three types of features were concatenated to form a 1×3072-dimensional composite feature vector. Only the target predicted values ​​(Kcat and Km, respectively) differed, and the training process remained consistent.

[0033] After model training, target enzyme sequences for optimization are selected, and functional sites are defined based on Swiss-Prot database annotations, Pfam domain analysis, and authoritative literature reports. A similarity threshold is assigned to each sequence, and the system searches for all protein sequences and their functional sites that meet the criteria. After the search, a combination of full-coverage random mutations and targeted mutations is used to generate mutant combinations: full-coverage random mutations perform 20 substitutions of amino acid residues at key sites, while targeted mutations perform residue substitutions within the user-defined functional site range. The instability index of each mutant sequence is calculated using the BioPython library; unstable sequences with an index ≥ 40 are removed, and candidate sequences with an index < 40 are retained.

[0034] Using the enhanced activity of two key enzymes in the vanillin synthesis pathway as an example, we conducted a virtual screening for highly efficient mutant enzymes: Feruloyl CoAsynthetase (Fcs) and Enoyl CoAhydratase / aldolase (Ech) play a catalytic role in the synthesis of vanillin, respectively, in the conversion of ferulic acid to feruloyl-CoA and its subsequent transformations. They are also involved in the generation of coenzyme A and acetyl-CoA. Therefore, improving the catalytic efficiency of these two enzymes is crucial to the overall yield.

[0035] A PubChem search revealed that FCS's UniProt ID is A0A2P1BT06, its substrate ferulic acid has a SMILES number of COC1=C(C=CC(=C1) / C=C / C(=O)O)O, and its product is feruloyl-CoA. CoA); Ech's UniProt ID is C3VA24, and its substrate feruloyl-CoA has the SMILES number CC1CCCN1C2=NC(=CC(=N2)N3CCN(CC3)C4=NC=C(C=C4C)C(=O)NCC5=CC=C(C=C5)F)C6=CC=C(C=C6)F, ultimately producing vanillin and acetyl-CoA. (CoA). After binding the Fcs sequence text to ferulic acid SMILES and the Ech sequence text to feruloyl-CoA SMILES, a similarity threshold of 0.35 was set in a Python script to query homologous sequences for both sequences. The results yielded 6 sequences similar to Fcs and 17 sequences similar to Ech. The collected similar sequence information includes binding regions and functional regions, identified by region labels and corresponding indices.

[0036] Then, a custom mutation method was defined, specifically setting the mutated residue index and target amino acid types. This example uses the system's default random mutation mode, performing random base substitutions on all functional sites. After mutation, a total of 3,000 potential Fcs optimized sequences and 7,317 potential Ech optimized sequences were generated. Their instability indices were calculated, and 807 sequences with unstable tendencies were removed from the Fcs potential optimized sequences, leaving 2,193 potential optimized sequences. After mutant construction, information such as mutation_id, mutation site combination, mutated residues, and instability index was uniformly recorded and saved in a CSV format for easy subsequent model retrieval.

[0037] Sequence processing module: Uses the ProtT5-xl model to extract semantic embeddings of protein sequences and uses average pooling to generate fixed vectors; Substrate processing module: Utilizes RDKit to generate Morgan fingerprints of substrate molecules and calls the MolT5 model to extract molecular semantic embeddings; Model call: All features are concatenated into a 3072-dimensional feature vector, which is then input into the ExtraTrees regression model to predict Kcat and Km respectively; Output: The generated Kcat / Km scores, InstabilityIndex, etc., will be used as subsequent filtering and sorting metrics.

[0038] Selection: Select sequences with InstabilityIndex < 40 and sort the results from highest to lowest by Kcat / Km.

[0039] The descriptions of the top 5 Fcs mutants by rating are shown below:

[0040] After removing 500 potentially unstable sequences from the Ech potential optimal sequences, 6817 potential optimal sequences remain. The top 5 are:

Claims

1. A method for predicting enzyme kinetic parameters and virtual screening based on multi-source feature fusion, characterized in that, Includes the following steps: Step 1: Obtain the enzyme sequence and its corresponding kinetic parameters Kcat and Km from public databases, perform data cleaning, redundancy removal and quality control, divide the training set and test set, and construct subsets of different intervals based on BLAST similarity; Step 2: Extract semantic embedding vectors of enzyme sequences using the ProtT5 protein pre-trained model; generate Morgan fingerprints of substrate molecules using the RDKit tool; extract semantic embedding vectors of substrates using the MolT5 molecular pre-trained model, and concatenate the obtained features to form a composite feature vector; Step 3: Train prediction models for Kcat and Km using extreme random tree regression models respectively, and improve model accuracy and generalization performance through cross-validation and hyperparameter tuning; Step 4: Based on the target enzyme sequence and functional site information, a combination of full-coverage random mutation and directed mutation is used to generate mutant combinations, calculate their instability index and remove unstable sequences. Step 5: Input the candidate mutants into the trained prediction model, predict the Kcat and Km values ​​in batches, calculate the Kcat / Km ratio, sort them according to the ratio, and output the candidate mutant sequences with the highest ratio.

2. The method as described in claim 1, characterized in that, In step 1, the training set and the test set are randomly shuffled proportionally. Based on the BLAST algorithm, a specified number of samples are randomly selected from the test set to construct four interval subsets with sequence similarity of 0%~40%, 40%~60%, 60%~80%, and 80%~100%, so that the data in each interval are evenly distributed in the test set.

3. The method as described in claim 1, characterized in that, The data cleaning in step 1 specifically includes: removing data with a Kcat value less than 10. - ³ s - ¹ or greater than 10 5 s - ¹, Km value less than 10 - ² μM or greater than 10 6 Abnormal data entries in μM; entries with invalid substrate SMILES encoding, enzyme sequence length less than 50 amino acid residues or greater than 1000 amino acid residues were removed; entries with missing key information, such as enzyme sequence, substrate SMILES, and kinetic parameters, were completely removed.

4. The method as described in claim 1, characterized in that, In step 2, the composite feature vector is a 1×3072-dimensional vector formed by concatenating the three types of features, which serves as the input to the prediction model.

5. The method as described in claim 1, characterized in that, In step 3, the ExtraTrees regressor in Scikit-learn is used as the prediction model. Hyperparameter tuning is achieved by combining 10-fold cross-validation and gridded hyperparameter search. The coefficient of determination R² is used as the evaluation index for model optimization.

6. The method as described in claim 1, characterized in that, In step 4, full-coverage random mutation involves 20 substitutions of amino acid residues at key sequence sites, forming a systematic mutation combination space; directed mutation involves residue substitution within the range of functional sites set by the user.

7. The method as described in claim 1, characterized in that, In step 4, the BioPython library is used to calculate the instability index of the mutant sequences. Mutant sequences with an instability index greater than or equal to 40 are removed, and mutant sequences with an instability index less than 40 are retained for subsequent screening.

8. The method as described in claim 1, characterized in that, The functional sites in step 4 are determined in the following ways: based on the manual annotation information of the Swiss-Prot database, the results of Pfam domain analysis, and the enzyme catalytic active site and substrate binding pocket related sites reported in authoritative literature, the range of functional sites is comprehensively defined; if the target enzyme has no known domain annotation, the conserved residue regions in the homologous sequence are identified as supplementary functional sites by using sequence conservation analysis tools.

9. The method as described in claim 1, characterized in that, In step 5, the sequences are sorted from largest to smallest according to the Kcat / Km ratio, and the top 10 candidate sequences with the highest ratios are selected as the final optimization results by default.

Images