A method for screening pancreatic lipase inhibiting peptides based on machine learning

Abstract

This invention discloses a method for screening pancreatic lipase inhibitory peptides based on machine learning, relating to the fields of bioinformatics and food science. The method includes the following steps: dataset construction and three-dimensional peptide structure characterization, statistical-based feature screening, outlier removal based on a dual-model consensus mechanism, feature recursive elimination and physicochemical significance optimization, data augmentation based on Gaussian noise injection, and final prediction model construction. This invention effectively solves the problems of insufficient feature representation, high noise interference, and weak model generalization ability in existing screening methods through multi-dimensional feature characterization, rigorous outlier cleaning, and data augmentation strategies for small samples. It can quickly and accurately screen highly active pancreatic lipase inhibitory peptides from massive peptide sequences, providing an efficient computational tool for the development of functional foods and lipid-lowering drugs.

Description

Technical Field

[0001] This invention relates to the fields of bioinformatics and food science, and in particular to a method for screening pancreatic lipase inhibitory peptides based on machine learning. Background Technology

[0002] Obesity and its associated metabolic syndromes (such as type 2 diabetes, hypertension, and cardiovascular disease) have become a global public health challenge. Pancreatic lipase (PL) is a key rate-limiting enzyme in the digestion of fat in the human intestine, responsible for hydrolyzing approximately 50%-70% of the triglycerides in the diet. By inhibiting the activity of pancreatic lipase, the absorption of dietary fat can be effectively reduced, thereby achieving the goal of weight control and obesity prevention.

[0003] Currently, orlistat, a pancreatic lipase inhibitor widely used in clinical practice, is effective, but long-term use often leads to gastrointestinal discomfort (such as steatorrhea and bloating), malabsorption of fat-soluble vitamins, and potential risks of liver and kidney damage. In contrast, food-derived bioactive peptides, due to their natural origin, high safety, easy absorption by the human body, and fewer side effects, have become a research hotspot for developing novel weight-loss functional factors.

[0004] However, existing technologies for screening and predicting pancreatic lipase inhibitory peptides still face the following significant technical bottlenecks and unresolved issues: (1) Traditional screening methods are inefficient: Traditional screening processes usually involve enzymatic digestion of proteins, complex separation and purification (such as ultrafiltration, gel chromatography, reversed-phase high-performance liquid chromatography), and cumbersome in vitro activity assays. This process is time-consuming, has high reagent costs, and often involves blind screening in a large number of unknown sequences, resulting in a low success rate and making it difficult to meet the needs of modern high-throughput screening.

[0005] (2) The experimental data is severely affected by noise: existing quantitative structure-activity relationship (QSAR) studies often directly cite IC from different literature. 50 Data sets are constructed using outliers. Due to differences in measurement conditions (such as substrate type, buffer pH, temperature, and enzyme concentration) between different laboratories, and the inevitable introduction of errors in experimental procedures, the original dataset often contains outlier samples. Most existing techniques lack rigorous outlier removal procedures and directly use noisy data for training, which seriously affects the model's generalization ability and prediction accuracy.

[0006] (3) Single feature representation dimension: Currently, most bioactive peptide prediction tools mainly rely on amino acid sequence features (such as amino acid composition, dipeptide composition, physicochemical properties, etc.). These features only treat peptides as a combination of characters, ignoring the three-dimensional spatial structure and topological properties of peptides as physical molecules in the microscopic environment. However, the binding of peptides to pancreatic lipase is essentially an intermolecular interaction, which depends on specific atomic arrangement, molecular volume, surface charge density, and physicochemical parameters such as van der Waals forces. It is difficult to capture these deep structure-activity relationships with only one-dimensional sequence information.

[0007] (4) High feature redundancy and lack of physical meaning: The number of peptide molecule descriptors is huge (up to thousands), containing a large number of redundant, highly collinear or activity-irrelevant features. Existing screening methods often fail to effectively remove highly correlated features, or simply rely on automatic selection by algorithms, ignoring the biophysical meaning behind the features. This results in the model performing reasonably well on the training set, but failing to truly reveal the key mechanisms of the "peptide-pancreatic lipase" interaction (such as spatial conformation, charge distribution, hydrogen bonding ability, etc.).

[0008] (5) Overfitting is prone to occur under small sample conditions: High-quality experimental data on pancreatic lipase inhibitory peptides are relatively scarce (usually only a few dozen to a few hundred). When the sample size is extremely small, the direct use of complex machine learning algorithms is prone to overfitting. Existing technologies lack effective data augmentation methods to expand the sample space, making it difficult to train robust and high-precision models.

[0009] In summary, there is an urgent need to develop a high-precision screening method that can effectively eliminate noisy samples, accurately capture key structure-activity characteristics (especially those based on three-dimensional structures), and overcome the limitations of small sample sizes through data augmentation techniques, in order to achieve efficient discovery of pancreatic lipase inhibitory peptides. Therefore, this invention proposes a machine learning-based method for screening pancreatic lipase inhibitory peptides to address the problems existing in the prior art. Summary of the Invention

[0010] To address the aforementioned problems, the present invention aims to propose a method for screening pancreatic lipase inhibitory peptides based on machine learning, thereby solving the problems of high data noise, high feature redundancy, and poor model overfitting and generalization ability caused by the scarcity of samples in the existing technology for screening pancreatic lipase inhibitory peptides.

[0011] To achieve the objectives of this invention, the present invention is implemented through the following technical solution: a method for screening pancreatic lipase inhibitory peptides based on machine learning, comprising the following steps: S1. Dataset Construction and Peptide 3D Structure Characterization Collect samples with clearly defined amino acid sequences and experimentally determined pancreatic lipase inhibitory activity values ​​(IC50). 50The initial dataset is constructed from peptide samples. Three-dimensional (3D) conformation files of each peptide sample are generated using peptide 3D structure prediction tools (such as PEP-FOLD). Then, molecular descriptor calculation software (such as PaDEL-Descriptor) is used to calculate molecular descriptors of all dimensions based on the 3D conformations, covering one-dimensional (1D), two-dimensional (2D) and three-dimensional (3D) features, forming the original feature matrix. S2. Initial feature screening based on statistics The original feature matrix obtained from S1 is preprocessed, and the specific steps are as follows: Missing value handling: Remove feature columns containing missing values; Variance screening: Calculate the variance of the remaining features and remove features with a variance of 0 (i.e., features whose values ​​do not change in all samples). Collinearity removal: Calculate the Pearson or Spearman correlation coefficient between features, set the correlation coefficient threshold to 0.9, and remove one of the feature pairs with a correlation coefficient greater than 0.9 to eliminate high redundancy between features and obtain the initial feature set; S3. Abnormal sample removal based on dual-model consensus mechanism To address potential noise interference in the experimental data, two heterogeneous algorithms, Random Forest (RF) and Extreme Gradient Boosting Tree (XGBoost, XGB), were employed for "double verification": Using initial feature set and IC 50 Activity value (or its negative logarithm pIC) 50 ), and establish RF regression model and XGB regression model respectively; K-fold cross-validation (e.g., 5-fold) was used to calculate the prediction residuals for each sample under the RF model and the XGB model respectively; Calculate the mean residual for each sample and set an outlier threshold (e.g., mean residual + 2 standard deviation). Samples with average residuals exceeding the anomaly threshold are marked as "abnormal samples" and removed. The remaining high-quality samples are retained as the cleaned training set. Anomalies are determined by the "consensus" of the two algorithms, avoiding accidental deletion due to the bias of a single model, and ensuring the objectivity and accuracy of data cleaning. S4. Feature recursion elimination and physicochemical significance optimization Deep feature engineering is performed on the cleaned training set to extract the key features that best characterize the inhibitory activity: Algorithm selection: Recursive feature elimination (RFECV) was performed based on RF and XGB algorithms respectively. The model performance of different feature subsets was evaluated by cross-validation, and two candidate feature sets were selected respectively. Manual selection: The two sets of candidate feature sets were merged and analyzed, and 14 core features that are highly related to the ligand-receptor interaction mechanism were manually screened in combination with their physicochemical significance (such as charge distribution, hydrogen bonding ability, steric hindrance, etc.).

[0012] The 14 core features specifically include: ATSC8i,SpMin3_Bhp,JGT,maxHBint6,AATSC3c,MATS3c,maxssCH2,ATSC3c,minHBint6,maxwHBa,MLFER_A,VE1_D,MDEC-33,ATSC4v. S5. Data Augmentation Based on Gaussian Noise Injection To address the issue of overfitting caused by a small sample size after cleaning (e.g., <100 samples), data augmentation is performed on the sample set obtained from S4: Keeping the target value (y) unchanged, random Gaussian noise is introduced into the feature matrix (X); The noise generation parameters are set as follows: mean is 0, and standard deviation is a certain proportion of the characteristic standard deviation (e.g., 0.05 times). The original samples are merged with noisy samples generated multiple times (e.g., 5 times) to construct the enhanced final training set. This step simulates the small perturbations in the experimental measurement, forcing the model to learn the essential rules of the features rather than memorizing specific values, which significantly improves the model's generalization ability. S6. Construction of the final prediction model Using the augmented training set obtained from S5, a random forest regression prediction model is constructed. The model is optimized by adjusting hyperparameters such as the number of decision trees (n_estimators), and the coefficient of determination (R²) of the model is evaluated using K-fold cross-validation. 2 The data augmentation operation, which involves the root mean square error (RMSE) and mean absolute error (MAE), is applied only to the split training set portion during each fold validation, while the validation set remains in its original state to ensure the objectivity of the evaluation results.

[0013] The beneficial effects of this invention are as follows: Extremely high data quality and model robustness: This invention innovatively proposes a dual-model consensus strategy of "random forest + extreme gradient boosting tree" for outlier removal. Compared with traditional 3σ principle or single model screening methods, this strategy utilizes two algorithms with different principles, Bagging and Boosting, for cross-validation. By calculating the average of the prediction residuals of the two models and setting an outlier threshold with the mean being twice the standard deviation, it can more accurately identify outliers caused by experimental errors and effectively avoid interference from noisy data on model training.

[0014] Overcoming the overfitting problem with small sample data: Addressing the technical challenge of scarce peptide activity experimental data (typically less than 100 samples), this invention introduces a data augmentation technique based on Gaussian noise. By keeping the target activity value constant, random noise with a mean of 0 and a standard deviation 0.05 times the feature standard deviation is applied to the feature matrix, and the original samples are merged with multiple generated noisy samples (in this embodiment, the sample size is expanded from 78 to 468), effectively expanding the sample space. The random forest model built on this basis achieves an average coefficient of determination R0 in five-fold cross-validation. 2 The model achieved a score of 0.9641 and maintained stable performance across all folds, demonstrating that it truly learned the intrinsic rules of structure-activity relationships rather than simply memorizing data, thus significantly improving the model's generalization ability.

[0015] The selected core features possess clear physical interpretability: the 14 core molecular descriptors ultimately selected in this invention are not random outputs of the algorithm, but rather have undergone dual verification through a combination of recursive feature elimination algorithm screening and artificial physicochemical meaning analysis. For example, maxHBint6 and minHBint6 characterize the intramolecular / intermolecular hydrogen bonding interaction capability of peptide molecules; the ATSC series descriptors characterize the autocorrelation between charge and spatial distribution; and descriptors such as JGT and VE1_D reflect the topological structure and geometric complexity of molecules. These features are directly related to the binding mechanism between peptides and pancreatic lipase active sites (such as electrostatic complementarity, hydrophobic interactions, hydrogen bond networks, etc.), providing a clear theoretical basis for subsequent molecular modification and de novo design of pancreatic lipase inhibitory peptides.

[0016] The prediction accuracy reaches industry-leading levels: the model constructed in this invention shows a high degree of agreement between the predicted values ​​and the actual experimental values ​​on the test set (R²≈0.96), and the prediction residuals follow a normal distribution and are randomly scattered around the zero axis, without exhibiting systematic bias. This prediction accuracy is significantly better than conventional quantitative structure-activity relationship (QSAR) models, enabling rapid and accurate screening of highly active pancreatic lipase inhibitory peptides from massive peptide sequences. It can effectively replace the cumbersome wet experimental screening process, significantly shorten the R&D cycle, reduce R&D costs, and provide an efficient computational tool for the development of functional foods and lipid-lowering drugs. Attached Figure Description

[0017] Figure 1 This is a schematic flowchart of the method for screening pancreatic lipase inhibitory peptides based on machine learning according to the present invention; Figure 2 This is the R-squared of the 5-fold cross-validation model in the embodiments of the present invention. 2 Score statistics chart; Figure 3 This is a model learning curve graph showing the increase in sample size in an embodiment of the present invention; Figure 4This is a linear fitting error diagram between the model prediction value and the actual experimental value in the embodiments of the present invention; Figure 5 This is a scatter plot of the residual distribution of the model training set and test set in an embodiment of the present invention; Figure 6 This is the QQ test plot of the normal distribution of the model prediction residuals in this embodiment of the invention; Figure 7 This is a beehive graph based on the importance and contribution of key features according to SHAP values ​​in this embodiment of the invention. Detailed Implementation

[0018] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0019] It should be noted that the technical means not described in detail in the following embodiments are all conventional means in the art, are not the key points of the invention, and will not be elaborated upon.

[0020] See Figure 1 This embodiment provides a method for screening pancreatic lipase inhibitory peptides based on machine learning, the method comprising the following steps: S1. Dataset Construction and Peptide 3D Structure Characterization The dataset in this embodiment was obtained from a systematic search of three major databases: Web of Science, Scopus, and CNKI, covering all years. To comprehensively obtain relevant literature, a search strategy combining Chinese and English keywords was adopted, with the main search terms including "peptide," "pancreatic lipase," "anti-obesity," and related terms.

[0021] The inclusion criteria require that the literature must meet the following conditions simultaneously: (1) Peptides with clear amino acid sequence information were reported; (2) The experimentally determined half-maximal inhibitory concentration (IC50) against pancreatic lipase is provided. 50 ) Specific numerical value.

[0022] The exclusion criteria removed the following types of data: (1) Samples containing non-standard amino acids or modified groups in the sequence to ensure consistency of molecular feature calculations; (2) IC 50Samples with inconsistent value units and lacking necessary information for effective conversion; (3) Literature whose experimental conditions (such as substrate type, reaction temperature, etc.) are not clearly described, which leads to doubts about the reliability of the data.

[0023] After the initial search, all obtained literature was manually checked and deduplicated, removing duplicate reports and data that did not meet the above criteria. Ultimately, 93 valid pancreatic lipase inhibitory peptide samples were collected and compiled, serving as the foundational dataset for subsequent model construction and analysis.

[0024] When building the dataset, IC 50 Value converted to negative logarithmic form pIC 50 As the target variable for model prediction, the feature matrix is ​​used in most existing peptide bioactivity prediction techniques, which tend to rely on macroscopic statistical features of peptide sequences, such as amino acid composition, dipeptide composition, or evolutionary information based on position-specific scoring matrices. While these traditional methods can capture sequence-level patterns, they often treat amino acids as single symbols or simple sets of physicochemical properties, making it difficult to deeply characterize the fine structural differences and microscopic electronic environment of peptide molecules at the atomic level. In contrast, this embodiment abandons traditional sequence statistical features and instead uses a more universal molecular descriptor system to characterize peptides. The innovation of this strategy is mainly reflected in the following three aspects: First, it achieves a leap from the "residue level" to the "atomic level." The molecular descriptors used in this embodiment are no longer limited to the statistical analysis of amino acid types, but directly quantify the microscopic physicochemical properties such as the topological distance between atoms within the molecule, charge distribution, ionization potential, and van der Waals volume. This allows for the keen detection of subtle activity changes caused by side chain modifications or stereoisomers, which is difficult to achieve with traditional sequence features.

[0025] Second, it breaks through the limitations of sequence length and homology. Traditional features based on sequence alignment or positional information are often limited by the uniformity of peptide length or depend on the abundance of known homologous sequences. In contrast, the universal molecular descriptor treats the peptide as a complete chemical entity, enabling unified and standardized mathematical characterization of peptides of arbitrary length and structure, which greatly expands the applicability and generalization ability of the model.

[0026] Third, it enhances the physical interpretability of structure-activity relationships. By introducing descriptors with clear physical meaning, the model can not only predict "how much activity there is," but also reveal "why it is active" from the perspectives of electrostatic interactions, hydrogen bond donor / acceptor capabilities, and steric hindrance, thereby revealing the deep molecular mechanism of the binding of pancreatic lipase inhibitory peptides to receptors.

[0027] In summary, this embodiment employs a universal molecular descriptor to construct a predictive model, which not only innovates existing feature extraction techniques but also provides a more precise microscopic perspective for the rational design and screening of highly active peptides. Specifically, the construction method involves using the PEP-FOLD3 server to predict the three-dimensional structure of the peptide and then using PaDEL-Descriptor software to calculate molecular descriptors based on this three-dimensional structure, obtaining a total of 1875 initial features encompassing 1D, 2D, and 3D dimensions, thus constructing the original feature matrix.

[0028] S2. Initial feature screening based on statistics To address the contradiction between the high-dimensional feature space resulting from general molecular descriptors and the relatively small number of peptide samples, this embodiment employs an unsupervised learning strategy for feature dimensionality reduction. First, using a variance thresholding method with a threshold of 0, descriptors with constant or minimally variable values ​​across all samples are eliminated, as these features do not carry information distinguishing activity differences. Subsequently, Spearman rank correlation coefficients are introduced to construct a feature correlation matrix. Considering the high degree of collinearity often present among peptide descriptors, which can lead to unstable model parameter estimation, this embodiment sets a correlation threshold of 0.90. For feature pairs with correlation coefficients higher than this threshold, only one is retained. This parameter setting is based on the principle of minimizing redundant features while retaining sufficient physicochemical information, reducing the computational complexity of subsequent modeling, and preventing overfitting, thus obtaining an initial feature set.

[0029] S3. Abnormal sample removal based on dual-model consensus mechanism To improve the quality of training data, this embodiment designs an outlier detection algorithm based on "dual-model consensus." It integrates Random Forest (RF) based on a Bagging strategy and Extreme Gradient Boosting (XGB) based on a Boosting strategy. In principle, RF improves stability by reducing variance, while XGB improves accuracy by reducing bias. The combination of the two can more comprehensively characterize the "normal" distribution of peptide structure-activity relationships. Regarding parameter settings, considering the small sample size, the RF model is set with 100 trees (n_estimators) and a maximum depth (max_depth) of 10; the XGB model is set with 100 trees, a maximum depth of 6, and a learning rate (learning_rate) of 0.1. Five-fold cross-validation (5-FoldCV) is used to obtain the predicted value of each sample under both models, and the absolute residual between the true value of each sample and the predicted values ​​of the two models is calculated. Finally, the average of the residuals of the two models is used as the judgment criterion, and the threshold is set to the mean + 2 times the standard deviation. Samples with residuals exceeding this threshold are considered "outliers," meaning their biological activity cannot be explained by current physicochemical laws, and are thus removed from the training set to ensure the robustness of the final model. This step removed 15 outlier samples, retaining 78 high-quality samples for subsequent feature engineering.

[0030] S4. Feature recursion elimination and physicochemical significance optimization In the feature selection stage, this embodiment abandons the one-sidedness of a single algorithm for screening and adopts a method combining recursive feature elimination and cross-validation (RFECV). This method iteratively builds the model, evaluates feature importance, and eliminates the weakest features until the optimal model performance is achieved. To ensure that the selected molecular descriptors have high universality and anti-interference ability, this embodiment runs two RFECV processes based on different principles: Random Forest and XGBoost. Random Forest focuses on the overall contribution of features, while XGBoost focuses on the ability of features to correct residuals. The final feature subset is the "intersection" of the two. The basis of this strategy is that only descriptors that are simultaneously identified as "key" by two completely different top algorithms are the core physicochemical factors that truly drive peptide activity. This greatly improves the reliability of the feature subset and avoids feature selection bias caused by algorithmic preferences. After feature selection and optimization, 14 key molecular descriptors were finally determined as input variables for the model: ATSC8i, SpMin3_Bhp, JGT, maxHBint6, AATSC3c, MATS3c, maxssCH2, ATSC3c, minHBint6, maxwHBa, MLFER_A, VE1_D, MDEC-33, and ATSC4v. These descriptors cover multiple dimensions such as autocorrelation, topological structure, atom type, and electrical properties. The specific calculation principles and physical meanings are explained below: To capture the interaction patterns between atoms within a molecule at specific distances, this embodiment selected several autocorrelation descriptors. Among them, ATSC8i (a Broto-Moreau autocorrelation descriptor with hysteresis 8th order centered based on ionization potential), ATSC3c (a Broto-Moreau autocorrelation descriptor with hysteresis 3rd order centered based on charge), and ATSC4v (a Broto-Moreau autocorrelation descriptor with hysteresis 4th order centered based on van der Waals volume) belong to the same family. They reflect the uniformity of the molecule's properties at different spatial scales by calculating the sum of the products of ionization potential, charge distribution, and van der Waals volume for atomic pairs separated by specific topological distances (8, 3, and 4 bond lengths, respectively). For example, ATSC3c can reveal the polarity pattern of local charge distribution in peptide chains, which is crucial for understanding the electrostatic interactions between peptides and enzyme active sites. Furthermore, AATSC3c (charge-based hysteresis 3rd order averaged centered Broto-Moreau autocorrelation descriptor) is an averaging process of ATSC3c, eliminating the influence of molecular size and focusing more on reflecting the intrinsic properties of charge distribution. MATS3c (charge-based hysteresis 3rd order Moran autocorrelation descriptor), on the other hand, introduces the Moran coefficient to measure the spatial autocorrelation between atomic charges separated by three bond lengths. A positive value indicates that closely spaced charges tend to cluster, while a negative value indicates an alternating distribution of charges. This is of great significance for analyzing the electrostatic potential surface of peptide molecules.

[0031] To describe the overall topological structure and connectivity of peptide molecules, descriptors such as SpMin3_Bhp and JGT are introduced. SpMin3_Bhp is a weighted descriptor based on Burden eigenvalues, specifically referring to the smallest eigenvalue (hysteresis 3) of the Burden matrix weighted by polarizability. The Burden matrix encodes atomic numbers, connectivity, and other information into its diagonal and off-diagonal elements. Its eigenvalues ​​can effectively characterize the branching of the molecule and the distribution of atomic polarizability in the topological structure, reflecting the induced dipole effect of the molecule in the microscopic environment. JGT (Global Topological Charge Index) is a comprehensive index that assesses the coupling degree between the overall charge distribution and topological structure of the molecule by calculating the topological distance along the charge transfer paths between all atomic pairs. It is often used to predict the global reactivity and stability of molecules. MDEC-33 (Modal Distance Counting Descriptor Based on Carbon Atoms) is specifically used to quantify the topological distance between two tertiary carbon atoms in a molecule. It reflects the compactness and steric hindrance effect of complex carbon structures in the peptide side chain or backbone, directly affecting the ease with which the peptide molecule enters the pancreatic lipase catalytic pocket. VE1_D (eigenvector coefficients based on the distance matrix) is derived from the distance matrix of the molecule. By calculating the sum of the coefficients of the eigenvectors of the distance matrix, the geometric complexity and spatial extensibility of the molecule as a whole are captured.

[0032] Finally, considering the crucial role of hydrogen bonding and specific atom types in receptor binding, this embodiment also selected descriptors related to hydrogen bonding ability and atom type. `maxHBint6` and `minHBint6` represent the strongest and weakest E-state descriptors that can form between atomic pairs separated by 6 topological distances within the molecule. These two indices directly quantify the potential of the peptide chain to form intramolecular hydrogen bonds during folding, which not only affects the secondary structure stability of the peptide but also determines its ability to form intermolecular hydrogen bonds with external receptors. `maxwHBa` (maximum weak hydrogen bond acceptor E-state index) focuses on the maximum electrical topological state value of the atom acting as a weak hydrogen bond acceptor in the molecule, typically involving atoms such as fluorine and sulfur or aromatic ring systems, which is indicative of identifying atypical hydrogen bond interaction sites. `maxssCH2` (maximum ssCH2 group E-state index) specifically characterizes the maximum electrical topological state of the methylene (-CH2-) group in the molecule, reflecting the electronic environment and reactivity of hydrophobic linker segments in the peptide chain. Furthermore, MLFER_A (linear free energy relation descriptor A) represents the overall hydrogen bond acidity of the molecule, i.e., the molecule's total ability to act as a hydrogen bond donor. This is a key physicochemical parameter determining whether a peptide molecule can form a stable bond with the nucleophilic group at the active site of pancreatic lipase. The comprehensive characterization using these 14 multi-dimensional descriptors lays a solid physicochemical foundation for the subsequent construction of a high-precision predictive model for pancreatic lipase inhibitory peptide activity.

[0033] S5. Data Augmentation Based on Gaussian Noise Injection The final prediction model was constructed using a random forest regressor, with the number of decision trees (n_estimators) set to 300 to ensure convergence and stability of the prediction results under the law of large numbers. Addressing the common "small sample size" problem in peptide drug development, this embodiment introduced Gaussian noise data augmentation technology. Specifically, the mean of the noise was set to 0, and the standard deviation was set to 5% of the standard deviation of each feature (i.e., a noise coefficient of 0.05). The augmentation factor was set to 5 times, generating 5 noisy samples for each original sample. The original data and the generated augmented data were merged to construct an augmented training set containing 468 samples. The principle behind this operation is to simulate small perturbations in the experimental environment, forcing the model to learn the essential patterns of data distribution rather than rote memorization of specific values, thereby significantly improving the model's generalization ability.

[0034] S6. Construction of the final prediction model A final random forest regression prediction model was constructed using the enhanced dataset. The number of trees in the model was set to 300, while other parameters remained at default or optimized states. The final model was evaluated using 5-fold cross-validation. The results showed that the model achieved a mean coefficient of determination (R²) of 0.9641 on the test set, with a significantly reduced root mean square error (RMSE), indicating that the model has extremely high prediction accuracy and generalization ability. Furthermore, the SHAP (SHapley Additive exPlanations) tool was used to calculate the SHAP values ​​of features and plot a beehive graph, visually demonstrating the positive or negative contribution of each feature to the prediction results, validating the biological rationale behind the model's decision-making. Through these steps, a complete construction of a high-precision activity prediction model from the original peptide sequence was achieved.

[0035] To verify the effectiveness of the prediction model based on dual-model consensus screening and data augmentation constructed in this embodiment, a multi-dimensional performance evaluation was performed on the trained model, as follows: 1. Coefficient of determination (R²) 2 ) Using the coefficient of determination R 2 R is a core metric for measuring the goodness of fit of a model. 2 This reflects the model's response to the dependent variable (peptide bioactivity pIC). 50 The variance explained by a model typically ranges from 0 to 1. A value closer to 1 indicates that the model's extraction of patterns from the molecular descriptors more accurately explains changes in activity.

[0036] 2. Root Mean Square Error (RMSE) To quantify the average deviation between predicted and actual values, the root mean square error (RMSE) is calculated. RMSE assigns a higher penalty weight to larger prediction errors, thus being more sensitive to outliers and rigorously reflecting the model's prediction accuracy. This metric is related to the target variable (pIC). 50 They have the same dimensions, making it easy to intuitively understand the physical magnitude of the error.

[0037] 3. Mean Absolute Error (MAE) As a supplement to RMSE, Mean Absolute Error (MAE) is also used for evaluation. MAE calculates the average level of the absolute value of the prediction error, reflecting the degree of deviation of the model's predicted values, and is unaffected by the direction (positive or negative) of the error. Compared to RMSE, MAE is less sensitive to outliers and more smoothly reflects the overall predictive ability of the model.

[0038] 4. Residual distribution and normality test Perform statistical tests on the model's prediction residuals. The residuals are defined as the difference between the actual and predicted values.

[0039] Residual distribution plot: Used to check whether the residuals are randomly distributed around the 0 axis. Ideally, the model residuals should show an irregular random distribution (homoscedasticity), indicating that the model has fully extracted the regularity in the data and has not missed any key systematic information.

[0040] QQ plot: Used to test whether the residuals follow a normal distribution. If the points in the plot approximately lie on a straight line, it indicates that the prediction error conforms to the Gaussian distribution assumption, verifying the statistical rationality of the model construction.

[0041] 5. SHAP value To address the "black box" problem of machine learning models, a game theory-based SHAP value is introduced as an evaluation metric for feature importance. The SHAP value quantifies the marginal contribution of each feature (molecular descriptor) to the prediction result for a specific sample. For a given feature, its SHAP value indicates the degree to which the presence of that feature causes the prediction result to deviate from the baseline value (average predicted value). A positive value indicates that the feature promotes increased peptide activity, while a negative value indicates inhibited activity. The core idea is to calculate the average marginal contribution of a feature across all possible feature combinations. Through this metric, this embodiment not only achieves accurate activity prediction but also reveals the intrinsic physicochemical mechanism of peptide structure-activity relationships.

[0042] In this embodiment, firstly, the coefficient of determination R² of the model is calculated using 5-fold cross-validation, such as... Figure 2 The cross-validation score plot shown indicates that the model's R-values ​​in five validations are... 2 The values ​​remained at a high level, with an average value stable at around 0.96, and the fluctuations in each validation result were small, demonstrating that the model has good robustness on different data subsets and that no obvious overfitting or underfitting was observed. Furthermore, as... Figure 3 The learning curves shown illustrate the changing trends of training set scores and cross-validation scores as the number of training samples increases. The two curves gradually converge and tend to be high, indicating that the aforementioned Gaussian noise data augmentation strategy effectively expands the sample space, enabling the model to fully learn data features and significantly improve generalization ability.

[0043] Next, the distribution of the model's prediction error was analyzed. For example... Figure 4 In the prediction error graph shown, the horizontal axis represents the experimentally determined true pIC. 50 The vertical axis represents the model's predicted values. The vast majority of sample points are closely distributed near the diagonal (i.e., the best-fit line), indicating a high degree of consistency between the predicted and actual values, and an extremely high linear fit. Meanwhile, as... Figure 5 The residual distribution plot shown indicates that the residuals of the training and test sets are randomly distributed around the zero axis, with most residual values ​​concentrated in a small range of ±0.5, exhibiting no obvious systematic bias. To further test whether the residuals conform to the normal distribution assumption, a plot is shown below. Figure 6 The residual QQ plot shown in the figure, with the sample points closely attached to the red reference line, confirms that the predicted residuals follow a normal distribution, further verifying the statistical reliability of the regression model.

[0044] Finally, to elucidate the biochemical mechanisms behind the model's predictions, the contribution of key features was visualized using the SHAP method. For example... Figure 7 The SHAP beehive diagram shown displays the selected key features arranged from top to bottom according to feature importance. Each point in the diagram represents a sample, and the color of the point represents the magnitude of the feature value (red for high values, blue for low values). The position of the point on the horizontal axis represents the impact of that feature on the model output (i.e., pIC). 50 The direction and extent of the influence of predicted values. The influence of molecular descriptors on the activity of pancreatic lipase inhibitory peptides can be divided into four categories: ① Autocorrelation descriptors related to charge distribution (ATSC8i, ATSC3c, AATSC3c, ATSC4v) ATSC8i, ranked first in importance, exhibits a clear positive correlation in its SHAP value distribution. ATSC8i characterizes the autocorrelation of ionization potentials among atomic pairs separated by eight chemical bonds in a molecule. A high ATSC8i value indicates that the peptide maintains a coordinated ionization potential distribution pattern over a long distance, which typically corresponds to the ordered folding of the peptide backbone and the regular arrangement of side chains. This structural feature facilitates precise spatial localization of the peptide near the active site of pancreatic lipase through long-range electrostatic interactions, thereby enhancing binding affinity. Conversely, a low ATSC8i value (blue dot) shows a negative contribution, suggesting that a disordered distribution of ionization potential may lead to a loose conformation of the peptide, making it difficult to form a stable enzyme-substrate complex. Similarly, ATSC3c and AATSC3c describe charge-based short-range autocorrelation. As can be seen from the SHAP plot, high values ​​for these two features also demonstrate a positive driving effect. This indicates that the charge aggregation or alternating distribution pattern in local regions of the peptide (within three bond lengths) is crucial to activity. For example, the orderly alternation of positive and negative charges may promote multi-point electrostatic complementarity between peptides and amino acid residues at the enzyme's active site (such as the Ser152, His263, and Asp176 catalytic triplet), thereby achieving efficient inhibition. Although ATSC4v is based on van der Waals volume calculations, its SHAP distribution also shows a high value corresponding to a positive contribution, indicating that moderate volume autocorrelation (i.e., coordinated distribution of side chain sizes) helps peptides embed into the enzyme's hydrophobic substrate channels.

[0045] ② Topological structure and molecular shape descriptors (SpMin3_Bhp, JGT, VE1_D, MDEC-33) SpMin3_Bhp (the minimum eigenvalue of the Burden matrix weighted by polarizability) exhibits a unique bidirectional influence pattern in the SHAP plot. Moderately high SpMin3_Bhp values ​​correspond to positive SHAP contributions, while extremely low or high values ​​produce negative effects. This inverted U-shaped relationship reveals a delicate balance in peptide polarizability distribution: moderate polarizability enhances induced dipole-dipole interactions, promoting peptide adsorption on the enzyme surface; however, excessively high polarizability may lead to overly flexible peptides, losing the rigid backbone required for entry into the active pocket. JGT (Global Topological Charge Index) quantifies the efficiency of charge transfer along topological paths. Figure 7 The results show that high JGT values ​​significantly contribute negatively to SHAP. This phenomenon can be explained by the fact that excessively high global charge transfer capacity may give the peptide excessive hydrophilicity or net charge, making it difficult for it to cross the hydrophobic "gated" region near the active site of pancreatic lipase, thus reducing the inhibition efficiency. VE1_D, as the sum of eigenvector coefficients based on the distance matrix, corresponds to a negative SHAP contribution with high values. This descriptor reflects the overall geometric complexity and spatial extensibility of the molecule. The negative correlation indicates that overly extended or complex peptide conformations may be unable to effectively enter the relatively narrow catalytic cavity of the enzyme due to steric hindrance. The SHAP distribution of MDEC-33 (tertiary carbon atom distance) also shows that high values ​​have a negative effect, confirming that excessive side chain branching or overly dense carbon backbone packing hinders the optimal geometric match between the peptide and the enzyme.

[0046] ③ Hydrogen bond interaction capability descriptors (maxHBint6, minHBint6, maxwHBa, MLFER_A) `maxHBint6` represents the strongest possible hydrogen bond interaction between atomic pairs separated by 6 topological distances within the molecule. The SHAP plot shows that a high value corresponds to a significant positive contribution, directly confirming the central role of hydrogen bonds in peptide-enzyme recognition. A strong internal hydrogen bond network not only stabilizes the active conformation of the peptide (such as β-turn or α-helical fragments) but also preserves sufficient hydrogen bond donor / acceptor sites for external interactions with the enzyme. Notably, `minHBint6` (the weakest hydrogen bond strength) shows a slight negative contribution trend, indicating that the presence of excessively weak hydrogen bond regions within the peptide may lead to conformational instability or exposure of unfavorable hydrophilic groups. The SHAP distribution of `maxwHBa` (the largest weak hydrogen bond acceptor E-state index) is relatively dispersed, but the overall trend shows that moderate weak hydrogen bond acceptor capacity (such as aromatic ring π electrons and sulfur atoms) can produce a positive effect. These non-classical hydrogen bonds may participate in CH-π or π-π stacking interactions with aromatic amino acids (such as Phe and Trp) or polar residues on the enzyme surface, enhancing binding specificity. A high MLFER_A (overall molecular hydrogen bond acidity) value corresponds to a negative SHAP contribution, suggesting that an excessively strong hydrogen bond donor capacity may lead to the formation of too many competitive hydrogen bonds between the peptide and the solvent water molecules, reducing its effective concentration at the enzyme active site.

[0047] ④ Specific atom type descriptors (maxssCH2, MATS3c) A high value of maxssCH2 (maximum methylene E-state index) indicates a positive SHAP contribution. The methylene group is the most common hydrophobic linker in polypeptide side chains, and its high E-state value signifies that the group is in a chemical environment with moderate electron density. This state facilitates close contact between the polypeptide and the hydrophobic substrate binding pocket of pancreatic lipase (composed of multiple hydrophobic residues such as Leu, Val, and Ile) via van der Waals forces, mimicking the fatty acid chain of the natural substrate triglyceride. MATS3c (charge-based Moran autocorrelation) describes the spatial autocorrelation of charges separated by three bond lengths. Its SHAP value distribution shows that positive correlation (charge aggregation) corresponds to a positive contribution, while negative correlation (charge alternation) contributes weakly. This suggests that local charge aggregation patterns (such as consecutive basic or acidic residues) may enhance Coulombic attraction with corresponding charge-complementary regions on the enzyme surface by forming "charge clusters."

[0048] The above SHAP value-guided mechanism analysis shows that highly efficient pancreatic lipase inhibitory peptides must simultaneously satisfy the following conditions: (1) The orderliness of long-range and short-range charge distribution (high ATSC series) ensures electrostatic complementarity; (2) Moderate topological complexity and polarization capability (moderate SpMin3_Bhp, low JGT) to ensure a balance between rigidity and flexibility in the conformation; (3) A strong internal hydrogen bond network (high maxHBint6) combined with appropriate external hydrogen bond donor and acceptor capabilities achieves dual optimization of conformational stability and ligand recognition; (4) Abundant hydrophobic linker units (high maxssCH2) promote tight binding with the enzyme's hydrophobic cavity. The revelation of these physicochemical laws provides a clear molecular modification strategy for the rational design of novel highly active inhibitory peptides.

[0049] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

Claims

1. A method for screening pancreatic lipase inhibitory peptides based on machine learning, characterized in that, Includes the following steps: S1. Collect peptide samples with clear amino acid sequences and experimentally determined pancreatic lipase inhibitory activity values, construct an initial dataset, use a peptide three-dimensional structure prediction tool to generate three-dimensional conformations of each peptide sample, and calculate molecular descriptors covering one-dimensional, two-dimensional and three-dimensional dimensions based on the three-dimensional conformations to form the original feature matrix. S2. Preprocess the original feature matrix formed in S1 to remove missing values, features with zero variance, and highly collinear features with correlation coefficients higher than a preset threshold, to obtain the initial screening feature set. S3. Using the initial screening feature set obtained in S2, establish a random forest regression model and an extreme gradient boosting tree regression model respectively. Calculate the prediction residuals of the samples under the two models through cross-validation, and identify and remove abnormal samples based on the consensus results of the two models to obtain the cleaned training set. S4. The recursive feature elimination algorithm is used to screen features in the cleaned training set in S3, and combined with molecular physicochemical analysis, the core features characterizing pancreatic lipase inhibitory activity are screened out. S5. Keep the target value of sample activity unchanged, inject Gaussian noise into the feature matrix corresponding to the core features selected in S4 to generate expanded samples, and merge the original samples and expanded samples to construct an enhanced training set. S6. Using the enhanced training set constructed in S5, a random forest regression prediction model is built for the screening and prediction of pancreatic lipase inhibitory peptides.

2. The method for screening pancreatic lipase inhibitory peptides based on machine learning according to claim 1, characterized in that: In step S1, the polypeptide three-dimensional structure prediction tool is the PEP-FOLD server, the molecular descriptor is calculated using PaDEL-Descriptor software, and the pancreatic lipase inhibitory activity value is IC50. 50 Value, converted to negative logarithm pIC 50 As the target variable for model prediction.

3. The method for screening pancreatic lipase inhibitory peptides based on machine learning according to claim 1, characterized in that: In step S2, the removal of highly collinear features uses Pearson or Spearman correlation coefficients, with a correlation coefficient threshold of 0.

9. For feature pairs with correlation coefficients higher than the threshold, one of the features is removed to eliminate feature redundancy.

4. The method for screening pancreatic lipase inhibitory peptides based on machine learning according to claim 1, characterized in that: In S3, the random forest regression model has 100 decision trees and a maximum depth of 10, the extreme gradient boosting tree regression model has 100 decision trees, a maximum depth of 6, and a learning rate of 0.1, and 5-fold cross-validation is used.

5. The method for screening pancreatic lipase inhibitory peptides based on machine learning according to claim 1, characterized in that: In S3, the specific process of removing abnormal samples is as follows: K-fold cross-validation is used to calculate the prediction residual of each sample under the random forest regression model and the extreme gradient boosting tree regression model respectively. The average value of the prediction residual of a single sample under the two models is calculated as the average residual. The abnormal threshold is set as the mean of the average residual of all samples plus twice the standard deviation. Samples with an average residual exceeding the abnormal threshold are marked as abnormal samples and removed.

6. The method for screening pancreatic lipase inhibitory peptides based on machine learning according to claim 1, characterized in that: In S4, there are a total of 14 core features, specifically: ATSC8i, SpMin3_Bhp, JGT, maxHBint6, AATSC3c, MATS3c, maxssCH2, ATSC3c, minHBint6, maxwHBa, MLFER_A, VE1_D, MDEC-33, and ATSC4v.

7. The method for screening pancreatic lipase inhibitory peptides based on machine learning according to claim 1, characterized in that: In S5, the parameters for Gaussian noise injection are: the noise mean is 0, and the standard deviation is 0.05 times the standard deviation of the corresponding feature. When performing cross-validation evaluation on the model, data augmentation based on Gaussian noise injection is applied only to the partitioned training set, while the validation set remains in its original state without augmentation.

8. The method for screening pancreatic lipase inhibitory peptides based on machine learning according to claim 1, characterized in that: In step S6, the number of decision trees in the random forest regression prediction model is set to 300, and the model is evaluated using 5-fold cross-validation. The evaluation metrics include the coefficient of determination R. 2 Root mean square error (RMSE) and mean absolute error (MAE).

Images