Corn microbial marker screening and phenotype prediction method

Abstract

The invention relates to a corn microbial marker screening and phenotype prediction method, and belongs to the technical field of corn breeding. The method comprises the following steps: preprocessing a microorganism sample: acquiring corn microorganism sample data, and performing data cleaning and standardization operation to generate a standardized microorganism abundance table; differential microorganisms are analyzed and screened through LEfSe, and finally a significant differential microorganism feature set is obtained; performing feature importance evaluation based on LightGBM: dividing the significant difference microorganism feature set screened by LEfSe into a training set and a test set; a LightGBM (LightGBM) model is constructed; the LightGBM model is used for calculating and evaluating the feature importance of the microorganisms; screening important characteristics to obtain a microbial marker for deeply analyzing key microbial community characteristics of corn phenotypic changes; and model verification and performance evaluation: verifying the prediction effectiveness of the screened microbial characteristics through model training, model verification and model performance evaluation. According to the method, the prediction effect of the model can be improved while important biological characteristics are reserved.

Description

Technical Field

[0001] The invention belongs to the technical field of corn breeding, and in particular relates to a method for screening corn microbial markers and predicting phenotypes. Background Art

[0002] In bioinformatics research, comparing different groups of biological samples to find significant differences is a key challenge. These differences may be due to experimental errors or real biological changes. It is crucial to screen out meaningful biological differences, especially the screening of microbial markers. Microbial markers can indicate specific environments, health states or diseases. By identifying and analyzing these markers, scientists can understand the health and disease states of the host, the impact of environmental changes on the microbial community, etc. However, the high-dimensional nature of microbiome data brings about problems such as data sparsity and multicollinearity. Feature selection helps to screen out useful information, reduce noise and redundancy, and improve the prediction ability of the model. The dynamic characteristics of the microbiome that change over time and environment have nonlinear effects on the host phenotype. Traditional methods find it difficult to capture these dynamic relationships, while complex machine learning methods can better reveal these patterns.

[0003] Linear discriminant analysis effect (LEfSe) is a common method for screening biomarkers in microbiome data. It measures the differences in microbial characteristics between groups through statistical tests (Kruskal-Wallis and Wilcoxon) and linear discriminant analysis (LDA). The results are easy to interpret, especially in biological studies. LEfSe can handle multi-group comparisons and evaluate the significant differences of microorganisms, but it assumes that the features are independent and it is difficult to capture complex nonlinear relationships and interactions between features, which limits its application in explaining the relationship between microbiome characteristics and phenotypic changes. In addition, LEfSe only considers statistical significance and ignores the contribution of features to model performance. LightGBM is a machine learning algorithm based on gradient boosting decision tree (GBDT), which is suitable for large-scale, high-dimensional data processing and is often used for classification, regression and feature selection tasks. Compared with LEfSe, LightGBM relies on model prediction performance, has poor interpretability, and is not easy to explain the relationship between features and target variables from a biological perspective. Although it can handle high-dimensional data, it may still be affected by noise and redundant features when faced with a large number of features, which increases model complexity and may lead to overfitting. Summary of the invention

[0004] In order to solve the above problems, the present invention provides a method for screening corn microbial markers and predicting phenotypes.

[0005] In order to achieve the above object, the present invention is implemented through the following technical solutions: A method for screening corn microbial markers and predicting phenotypes comprises the following steps: S1. Microbial sample pretreatment: Obtain corn microbial sample data, perform data cleaning and standardization operations, generate a standardized microbial abundance table, perform Alpha diversity analysis to understand the internal structure of the microbial community, and verify whether the pretreatment method is reasonable; S2.LEfSe analysis to screen differential microorganisms: group the samples in the microbial abundance table; use the Kruskal-Wallis rank sum test to screen the microbial features with significant differences; use the Wilcoxon rank sum test to compare the two groups to screen the microbial features with stable and significant differences between the groups; use linear discriminant analysis to calculate the effect size of microbial features with stable and significant differences, and finally obtain the set of microbial features with significant differences; S3. Feature importance evaluation based on LightGBM: Divide the significantly different microbial feature set screened by LEfSe into a training set and a test set; construct a LightGBM model; the LightGBM model performs computational evaluation on the importance of microbial features; screen important features to obtain microbial markers for in-depth analysis of key microbial community characteristics of maize phenotypic changes; S4. Model validation and performance evaluation: The predictive validity of the screened microbial features was verified through model training, model validation and model performance evaluation.

[0006] Furthermore, step S1 analyzes the microbial community in the corn microbial sample data by high-throughput sequencing technology to obtain a microbial abundance table, wherein the high-throughput sequencing technology includes 16S rRNA gene sequencing and metagenomic sequencing.

[0007] Further, step S2 includes the following steps: S21. Input data and grouping: Group the corn samples in the microbial abundance table according to the experimental design. Each sample has a corresponding classification label to distinguish the differences between the groups and ensure that the subsequent analysis can accurately capture the differences in microbial characteristics of different groups; the classification label is the corn phenotype, which includes genotype and growth stage; S22. Kruskal-Wallis rank sum test to screen significant microbial features: LEfSe first performed a Kruskal-Wallis rank sum test on the microbial abundance table; the significance level was set to a p value less than 0.05, and microbial features with statistically significant differences between different groups were screened out. The formula is as follows: , in, H Represents the test statistic, which is used to determine whether there is a significant difference between groups. N represents the total number of samples, m Indicates the number of groups, ni represents the number of samples in the i-th group, R i represents the rank sum of the i-th group; S23. Wilcoxon rank sum test was used to compare the two groups to screen out the microbial characteristics with stable and significant differences between the groups: For the microbial characteristics with significant differences through Kruskal-Wallis test, Wilcoxon rank sum test was further performed to evaluate the specific differences between each group; S24. Linear discriminant analysis was used to calculate the effect size of microbial features with stable and significant differences: Linear discriminant analysis (LDA) was used to calculate the effect size of microbial features with stable and significant differences. By setting the effect size threshold, the microbial features that contributed more to the differences between groups were screened out, and then a set of microbial features with significant differences was obtained. The effect size calculation formula of LDA is: , Among them, LDA Score represents the effect size of microbial characteristics. The larger the value, the greater the contribution of this microorganism to the difference. μ g1 , μ g2 They represent the characteristic means of the first group and the second group respectively, and σ represents the overall standard deviation of the characteristic.

[0008] Further, step S3 includes the following steps: S31. Data preparation: The significantly different microbial feature set screened by LEfSe is used as the input feature of the LightGBM model, and the corn phenotype is used as the prediction target of the model; the significantly different microbial feature set screened by LEfSe is divided into a training set and a test set; S32. Construct LightGBM model: Use the histogram-based decision tree construction method to discretize the abundance values ​​of microorganisms in the microbial abundance table into several intervals and construct a histogram. The model splits by selecting the feature that produces the largest reduction in the loss function, generates binary decision trees layer by layer, and finally obtains the LightGBM model. Set the basic parameters of the LightGBM model, including the number of leaf nodes num_leaves, learning rate learning_rate, and maximum depth max_depth of the tree for each tree; S33. Feature importance evaluation: The LightGBM model provides a variety of feature importance evaluation methods, including feature importance calculation based on gain, coverage, and split times. The significant difference microbial feature set is evaluated through feature importance evaluation to obtain the importance score of the significant difference microbial feature. S34. Screening important features: sorting the significantly different microbial features according to the importance scores of the significantly different microbial features obtained in step S33, and screening important features as microbial markers; S35. Predict corn phenotype based on the microbial markers obtained in step S34 to obtain the final predicted phenotype.

[0009] Further, step S4 includes the following steps: S41: Model training: Use the significantly different microbial features in the training set as the input features of the LightGBM model, take the corn phenotype as the prediction target of the model, train the LightGBM model on the training set, learn the association pattern between the significantly different microbial features and phenotypic changes, and finally obtain a trained LightGBM model; S42: Model validation: Validate the trained LightGBM model on the test set to evaluate the prediction performance of the model on unseen data; obtain the overall performance indicators of the model through the prediction results of the test set; use the cross-validation method to further verify the robustness of the model; S43: Model performance evaluation: Combine the test set and cross-validation results to comprehensively evaluate the performance indicators of the model.

[0010] The advantages of the present invention are: The present invention uses LEfSe statistical tests and LDA to screen microbial features to ensure that the screened features have significant differences between different groups. This screening method has strong biological interpretability and is consistent with the understanding of the data by experts in the field. After LEfSe dimensionality reduction, the feature set trained by the LightGBM model is more interpretable, and the analysis results are more consistent with biological and statistical significance. This helps to more easily explain the contribution of microbial features to the phenotype in model applications, and improves the practicality and convincingness of the model. Preliminary dimensionality reduction of microbiome data by LEfSe can significantly simplify the data feature space, improve the training efficiency and predictive performance of the LightGBM model, and ensure the interpretability of the model. This combination method not only reduces the complexity of the model in the analysis of high-dimensional microbiome data, but also retains the biomarkers that have the greatest impact on the target phenotypic changes.

[0011] Many features in the microbiome data may have a weak or no correlation with the target variable (such as host genotype), or may even be noisy. Although LightGBM itself can handle a large number of features, redundant or irrelevant features will increase the computational cost and may affect the performance of the model. LEfSe can effectively remove these redundant features, helping LightGBM focus on the most important features, thereby improving the model's prediction accuracy. The present invention combines LEfSe and LightGBM. LEfSe first reduces the dimensionality to screen important features, and then LightGBM models based on these features. It can improve the prediction effect of the model while retaining important biological features. BRIEF DESCRIPTION OF THE DRAWINGS

[0012] The accompanying drawings are used to provide further understanding of the present invention and constitute a part of the specification. They are used to explain the present invention together with the embodiments of the present invention and do not constitute a limitation of the present invention.

[0013] Figure 1 is a flow chart of the method of the present invention; Figure 2 This is the changing trend of the cross-validation error of the genotype binary classification model; Figure 3 This is a graph showing the relationship between the number of features and the cross-validation loss; Figure 4 This is the confusion matrix diagram for maize genotype phenotype prediction based on LightGBM; Figure 5 This is the ROC curve diagram for maize genotype phenotype prediction based on LightGBM; Figure 6 This is the changing trend of the cross-validation error of the multi-classification model in the growth stage; Figure 7 This is the relationship between the number of model features and cross-validation loss in the growth phase; Figure 8 This is the confusion matrix diagram for predicting corn growth stage phenotype based on LightGBM; Fig. 9 This is the ROC curve for predicting corn growth stage phenotype based on LightGBM. DETAILED DESCRIPTION

[0014] The following will be combined with the drawings in the embodiments of the present invention to clearly and completely describe the technical solutions in the embodiments of the present invention. Obviously, the described embodiments are only 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.

[0015] Example 1 In this embodiment, Figure 1 As shown, this embodiment provides a method for screening corn microbial markers and predicting phenotypes, and the specific steps include: S1. Microbial sample pretreatment: Obtain corn microbial sample data, perform data cleaning and standardization operations, generate a standardized microbial abundance table, perform Alpha diversity analysis to understand the internal structure of the microbial community, and verify whether the pretreatment method is reasonable.

[0016] Specifically, the microbial community in the corn microbial sample data is analyzed by high-throughput sequencing technology to obtain a microbial abundance table, and the high-throughput sequencing technology includes 16S rRNA gene sequencing and metagenomic sequencing.

[0017] S2.LEfSe analysis to screen differential microorganisms: group samples in the microbial abundance table; use the Kruskal-Wallis rank sum test to screen microbial features with significant differences; use the Wilcoxon rank sum test to compare two groups to screen out microbial features with stable and significant differences between groups; use linear discriminant analysis to calculate the effect size of microbial features with stable and significant differences, and finally obtain a set of significantly different microbial features.

[0018] Specifically, step S2 includes the following steps: S21. Input data and grouping: Group the corn samples in the microbial abundance table according to the experimental design. Each sample has a corresponding classification label to distinguish the differences between the groups and ensure that the subsequent analysis can accurately capture the differences in microbial characteristics of different groups; the classification label is the corn phenotype, which includes genotype and growth stage; S22. Kruskal-Wallis rank sum test to screen significant microbial features: LEfSe first performed a Kruskal-Wallis rank sum test on the microbial abundance table; the significance level was set to a p value less than 0.05, and microbial features with statistically significant differences between different groups were screened out. The formula is as follows: , in, H Represents the test statistic, which is used to determine whether there is a significant difference between groups. N represents the total number of samples, m Indicates the number of groups, n i represents the number of samples in the i-th group, R i represents the rank sum of the i-th group; S23. Wilcoxon rank sum test was used to compare the two groups to screen out the microbial characteristics with stable and significant differences between the groups: For the microbial characteristics with significant differences through Kruskal-Wallis test, Wilcoxon rank sum test was further performed to evaluate the specific differences between each group; S24. Linear discriminant analysis was used to calculate the effect size of microbial features with stable and significant differences: Linear discriminant analysis (LDA) was used to calculate the effect size of microbial features with stable and significant differences. By setting the effect size threshold, the microbial features that contributed more to the differences between groups were screened out, and then a set of microbial features with significant differences was obtained. The effect size calculation formula of LDA is: , Among them, LDA Score represents the effect size of microbial characteristics. The larger the value, the greater the contribution of this microorganism to the difference. μ g1 , μ g2 They represent the characteristic means of the first group and the second group respectively, and σ represents the overall standard deviation of the characteristic.

[0019] S3. Feature importance evaluation based on LightGBM: Divide the significantly different microbial feature set screened by LEfSe into a training set and a test set; construct a LightGBM model; use the LightGBM model to perform computational evaluation on the importance of microbial features; screen important features to obtain microbial markers for in-depth analysis of key microbial community characteristics of maize phenotypic changes.

[0020] Specifically, step S3 includes the following steps: S31. Data preparation: The significantly different microbial feature set screened by LEfSe is used as the input feature of the LightGBM model, and the corn phenotype is used as the prediction target of the model; the significantly different microbial feature set screened by LEfSe is divided into a training set and a test set; S32. Construct LightGBM model: Use the histogram-based decision tree construction method to discretize the abundance values ​​of microorganisms in the microbial abundance table into several intervals and construct a histogram. The model splits by selecting the feature that produces the largest reduction in the loss function, generates binary decision trees layer by layer, and finally obtains the LightGBM model. Set the basic parameters of the LightGBM model, including the number of leaf nodes num_leaves, learning rate learning_rate, and maximum depth max_depth of the tree for each tree; S33. Feature importance evaluation: The LightGBM model provides a variety of feature importance evaluation methods, including feature importance calculation based on gain, coverage, and split times. The significant difference microbial feature set is evaluated through feature importance evaluation to obtain the importance score of the significant difference microbial feature. S34. Screening important features: sorting the significantly different microbial features according to the importance scores of the significantly different microbial features obtained in step S33, and screening important features as microbial markers; S35. Predict corn phenotype based on the microbial markers obtained in step S34 to obtain the final predicted phenotype.

[0021] S4. Model validation and performance evaluation: The predictive validity of the screened microbial features was verified through model training, model validation and model performance evaluation.

[0022] Specifically, step S4 includes the following steps: S41: Model training: Use the significantly different microbial features in the training set as the input features of the LightGBM model, take the corn phenotype as the prediction target of the model, train the LightGBM model on the training set, learn the association pattern between the significantly different microbial features and phenotypic changes, and finally obtain a trained LightGBM model; S42: Model validation: Validate the trained LightGBM model on the test set to evaluate the prediction performance of the model on unseen data; obtain the overall performance indicators of the model through the prediction results of the test set; use the cross-validation method to further verify the robustness of the model; S43: Model performance evaluation: Combine the test set and cross-validation results to comprehensively evaluate the performance indicators of the model.

[0023] Example 2 In this example, in order to explore the differences in the temporal dimension and genotype between the rhizosphere microbial communities of B73 and Mo17, an experiment was conducted in three adjacent fields in Xinxiang City, Henan Province.

[0024] 1. Data preprocessing In each field, B73 and Mo17 maize varieties were planted alternately, and rhizosphere soil samples were systematically collected at four key growth stages: week 4 (DAS28), week 6 (DAS42), week 8 (DAS56), and week 10 (DAS70). In this study, a total of 72 samples were collected, covering two maize varieties, three experimental plots, and four growth periods, with three replicates in each period.

[0025] All samples were subjected to metagenomic sequencing, and after strict quality control, 70 high-quality rhizosphere samples were obtained, containing a total of 2,474,420,040 sequencing reads for detailed microbial community analysis; after taxonomic annotation of the sequencing results, a microbial abundance table was obtained, in which bacteria were absolutely dominant in the maize rhizosphere microbiome, accounting for more than 91% of all taxonomic units. These bacteria belonged to 49 phyla, 104 classes, 231 orders, 539 families, 2025 genera, and 6381 species.

[0026] In order to eliminate the possible bias introduced by the difference in sequencing depth (i.e., the total number of reads obtained from sequencing each sample) of different samples, this example dilutes the sample data. The purpose of dilution is to standardize the sequencing depth of all samples to the same level to ensure the comparability of data between samples. In this example, the sequencing depth of all samples is uniformly diluted to 100,000 reads to achieve standardization of sequencing depth and provide unified basic data for subsequent analysis.

[0027] The data presented here were at the genus level, and 35 samples of each genotype (B73 and Mo17) were analyzed at four growth stages: DAS28 (n=8), DAS42 (n=9), DAS56 (n=9), and DAS70 (n=9).

[0028] The present invention compared the Alpha diversity of B73 and Mo17 and found that the diversity of rhizosphere microbial communities of B73 and Mo17 showed temporal differences at different growth stages with observable changes, which is a prerequisite for subsequent feature selection and machine learning modeling.

[0029] Table 1 Comparative analysis of Alpha diversity index between B73 and Mo17 As shown in Table 1, Alpha diversity refers to the diversity of microbial communities in a single sample, usually measured by indicators such as species richness (number) and uniformity (distribution). Shannon index: measures the diversity of a community, taking into account species richness and uniformity. The higher the value, the more complex and evenly distributed the community is. Simpson index: measures the uniformity of a community, reflecting the probability that a certain species dominates the community. The smaller the value, the lower the monopoly of the dominant species. Pielou index: reflects the uniformity of species distribution in a community. Relative to the maximum possible diversity, the closer the value is to 1, the more uniform the species distribution is. The values ​​in the table represent the mean ± standard deviation (SD) of the diversity index. The P value is used to indicate the possibility that a certain result is due to chance. The smaller the value, the less likely it is that the result is due to chance, but the greater the possibility that there is a real difference. The significance of the difference between groups was calculated by the Wilcoxon rank sum test and expressed as a P value.

[0030] For the Alpha diversity index, the significance of the differences between the groups was determined by the Wilcoxon rank sum test to see the differences between the two genotypes, and the values ​​represent the mean ± standard deviation (SD) of the diversity index.

[0031] Table 2 Comparative analysis of Alpha diversity index between B73 and Mo17 at four growth stages For the Alpha diversity index, the significance of the differences between groups was determined by Dunn's Kruskal-Wallis multiple comparison test, which was used to look at the differences between multiple groups (because we had four growth stages), as shown in Table 2. In the results, the diversity index is expressed as mean ± standard deviation (SD), and the significant differences between groups are annotated by letter marks (such as a, b, c, d): different letters indicate significant differences between groups, and the same letters indicate insignificant differences. DAS28 represents the fourth week, DAS42 represents the sixth week, DAS56 represents the eighth week, and DAS70 represents the tenth week.

[0032] 2. The following is to screen biomarkers and predict phenotypes from the perspective of genotype and growth stage. (1) LEfSe analysis and screening of differentially expressed microorganisms (genotypes) To further understand the differences in the microbial communities between B73 and Mo17, we performed LEfSe analysis. LEfSe performs Kruskal-Wallis and Wilcoxon tests in-house, simplifying the complex statistical process into an LDA score. Features with high LDA scores are not only significantly different, but also have large effect sizes, thereby accurately identifying features with significant statistical differences in rhizosphere microorganisms. There are 35 samples each for the two genotypes B73 and Mo17.

[0033] The analysis results showed that there were 89 genera with LDA scores greater than 2, 32 genera with LDA scores greater than 2.5, and 12 genera with LDA scores greater than 3, as shown in Table 3.

[0034] Table 3 LDA score table of rhizospheric bacterial genera of B73 and Mo17 (LDA>3) As shown in Table 3, this table is used to show the significant differences of different bacterial genera in the rhizosphere bacteria of B73 and Mo17 varieties, as well as the relative importance of these differences (LDA score>3). The higher the LDA score, the greater the contribution of the bacterial genus in distinguishing the rhizosphere microbial communities of B73 and Mo17.

[0035] (2) Feature importance evaluation based on LightGBM (genotype) The 89 genera with LDA scores greater than 2 are converted into the format of LightGBM, and then the parameters are set and 10-fold cross validation is performed. This cross-validation method randomly divides the dataset into 10 subsets and iterates the training 10 times, with 9 copies used for training and 1 copy used for validation each time to ensure the robustness of the model under different data partitions. Through cross-validation, LightGBM calculates the loss value on the test set after each round of iteration and determines the optimal number of iterations to avoid overfitting. After determining the optimal number of iterations, the final model is trained under the optimal iteration using the full data. Next, the feature importance is extracted and a list of features sorted by importance is returned. Feature importance is presented in three indicators: information gain (Gain), coverage (Cover) and split frequency (Frequency). The feature importance table of the genotype binary classification model is shown in Table 4. Figure 2This is the trend of the cross-validation error of the genotype binary classification model. As can be seen from the figure, the curve gradually decreases with the increase in the number of iterations from the initial iteration, indicating that the model significantly reduced the error by adjusting the weights in the early stage of training. In the later stage, the curve tends to be stable, and may even rise slightly, which indicates that the model may have reached the best iteration round, and the subsequent training effect has limited improvement or even begins to overfit. The main purpose of this figure is to help determine the best training round to balance the performance and training efficiency of the model. At the same time, it can also intuitively display the learning status of the model, so as to determine whether there is underfitting or overfitting.

[0036] Table 4. Feature importance table of genotype binary classification model (top 20%) Next, a feature screening experiment was conducted to find out the most important number of features for model prediction by gradually increasing the number of features used in the model and observing its effect on the cross-validation loss (Log Loss). The whole process was repeated 10 times. In each experiment, the number of features was gradually increased from the most important one until all features were included. Each experiment used leave-one-out cross-validation (LOOCV), that is, one sample was left out from all samples as a test set each time, and the remaining samples were used for training, thereby ensuring the robustness of the model under different data combinations.

[0037] The purpose of repeating the experiment 10 times is to reduce the random impact of data partitioning in a single experiment. Through repeated cross-validation, we can observe the overall trend of the change in the number of features on the model loss and obtain more robust evaluation results. Finally, we average the results of the 10 experiments to eliminate the interference of accidental factors, making the trend of the impact of the number of features on model performance clearer, thereby helping to identify the optimal number of features with the lowest loss in multiple experiments.

[0038] Figure 3 The figure shows the relationship between the number of features and the cross-validation loss. It can be seen from the figure that when 6 features are used, the cross-validation error curve stabilizes. Therefore, we define these 6 genera as microbial markers. These 6 genera are Sphingomonas, Fibrobacter, Xinella, Sandobacillus, Lantella and Gemella; (3) Phenotype prediction - model validation and performance evaluation (genotype) Using the 6 selected microbial markers as input features (independent variables) and the corn variety corresponding to each sample as the prediction target (dependent variable), we built a classification model based on LightGBM to predict the corn genotype of the sample. A loop was established, in which one sample was set aside as a test set each time, and the remaining samples were used to train the model. In each loop, after removing the current sample from the dataset, a 10-fold cross-validation was performed with the remaining dataset, and the best number of iterations was recorded. The final model was then trained at this number of iterations, and the trained model was used to predict the remaining test samples. After the entire loop is completed, the prediction results of each sample are saved, including whether the prediction is correct, the actual label, and the predicted label.

[0039] Figure 4 This is a confusion matrix diagram for maize genotype phenotype prediction based on LightGBM, which helps us intuitively understand the classification effect of the model. The count of each cell in the confusion matrix shows the number of times a specific actual category is predicted as a certain category. The higher the value on the diagonal, the better the classification effect of the model. The values ​​outside the diagonal represent the number of classification errors, and the lower the value, the better. By looking at the non-diagonal part, you can understand the tendency of the model to misclassify some actual categories as other categories. The overall classification effect of the model is good, and Mo17's recognition is slightly better than B73, but the difference is not large.

[0040] In order to fully understand the classification performance based on the LightGBM model, we plotted the ROC curve and listed the main classification performance indicators, such as Figure 5 As shown. Based on the prediction probability output by the model, different TPRs and FPRs were obtained by gradually adjusting the prediction threshold, thereby generating the ROC curve. The AUC value is close to 1, indicating that the model has a high discrimination ability in distinguishing B73 and Mo17 varieties. The large AUC area indicates that the LightGBM model can effectively distinguish the genotype differences of corn and has a high prediction accuracy.

[0041] To more comprehensively evaluate the model performance, we calculated common classification metrics and compared them with models under different feature combinations, including the following metrics: Accuracy: The overall prediction accuracy of the model for all samples. Precision: The proportion of samples correctly predicted as positive by the model to all samples predicted as positive.

[0042] Recall: The proportion of samples that the model correctly predicts to be positive to all samples that are actually positive. (Accuracy of predicting B73) F1 Score: The harmonic average of precision and recall, used to measure the overall performance of the model. Specificity: The proportion of negative classes that the model correctly identifies. (Predicting the accuracy of Mo17) Kappa coefficient: measures the gap between model classification effect and random classification, reflecting consistency.

[0043] Table 5 Model evaluation table (genotype) As shown in Table 5, "all genera" includes 2025 genera obtained by annotation; "P value genera" are the first 89 genera selected according to the P value ranking of the Wilcoxon rank sum test; "random genera" are 89 randomly selected genera; "LDA>2" is 89 differential microbial genera with LDA scores greater than 2 screened out by LEfSe; "LDA>2.5" is 32 genera with LDA scores greater than 2.5 screened out by LEfSe; "LDA>3" is 12 genera with LDA scores greater than 3 screened out by LEfSe; the final "microbial markers" are 6 key genera obtained based on LightGBM feature screening. The model of the microbial marker feature set performed well in all indicators and was significantly better than other feature selection methods. This shows that the 6 key genera screened based on LightGBM can effectively distinguish B73 and Mo17 and have high classification performance. In contrast, the model containing all 2025 genera and the model of 89 randomly selected genera had lower accuracy and Kappa coefficient (0.800 and 0.600), indicating that the unscreened feature set was less effective in accurately distinguishing the two varieties. Among other feature selection methods, 32 genera with LDA scores greater than 2.5 performed relatively well in terms of recall (0.853) and specificity (0.857), but were still lower than the results of the "microbial marker" feature set, indicating that the use of a small number of important features after screening can improve the specificity of classification. In addition, we found that Mo17 was easier to predict, reflecting that the microbial community structure of Mo17 at different growth stages may be more stable or consistent, with less variability.

[0044] (4) LEfSe analysis and screening of differential microorganisms (growth stage) We performed LEfSe analysis at the genus level for all samples in four growth stages to identify statistically different microorganisms. The fourth week (DAS28) contained 16 samples, the sixth week (DAS42) contained 18 samples, the eighth week (DAS56) contained 18 samples, and the tenth week (DAS70) contained 18 samples.

[0045] The analysis results show that there are 75 genera with LDA scores greater than 2, 38 genera with LDA scores greater than 2.5, and 17 genera with LDA scores greater than 3 (see Table 6).

[0046] Table 6 LDA scores of rhizospheric bacteria in B73 and Mo17 at four stages (LDA>3) DAS28 represents the fourth week, DAS42 represents the sixth week, DAS56 represents the eighth week, and DAS70 represents the tenth week.

[0047] (5) Feature importance evaluation based on LightGBM (growth stage) After converting the 75 genera with LDA scores greater than 2 into LightGBM format, set the parameters and perform 10-fold cross validation to ensure the robustness of the model under different data partitions. Through cross validation, LightGBM calculates the loss value on the test set and determines the optimal number of iterations to avoid overfitting. After determining the optimal number of iterations, the final model is trained under this iteration using the full dataset, and the feature importance is extracted to generate a list of features sorted by importance. Feature importance is presented through three indicators: information gain (Gain), coverage (Cover), and split frequency (Frequency). See Table 7.

[0048] Table 7. Importance of multi-classification model features in the growth stage (top 20%) Figure 6 This is the trend of the cross-validation error of the multi-classification model in the growth stage. It can be seen from the figure that the curve gradually decreases with the increase in the number of iterations in the initial stage, indicating that the model significantly reduces the error by adjusting the weights in the early stage of training. However, after reaching the best iteration round, the curve begins to rise gradually, indicating that the model has overfitted when continuing training. Overfitting will cause the performance of the model on the test set to deteriorate, indicating that it overfits the training data. However, this does not affect the final model performance, because we choose the best iteration round as the training result of the model, which effectively avoids the negative impact of overfitting.

[0049] Next, a feature screening experiment was conducted to gradually increase the number of features used in the model and observe its impact on the cross-validation loss (Log Loss) to identify the number of features that are most important for model prediction. Each experiment started with the most important features and gradually increased them, using leave-one-out cross-validation (LOOCV) to ensure the robustness of the model. To reduce the impact of randomness, the experiment was repeated 10 times and the average was taken to obtain the overall trend of the change in the number of features on model performance and determine the optimal number of features with the lowest loss.

[0050] like Figure 7As shown in the figure, when 11 features are used, the cross-validation error curve tends to be stable, so these 11 genera are defined as microbial markers, namely, Rhodothricinus, Cellulobacterium, Sandmonas, Sinorhizobium, Near-acidobacterium, Acidobacterium, Box-shaped Coccus, Grass Spirillum, Leaf Bacillus, Verrucomicrobia and Flavobacillus.

[0051] (6) Phenotype prediction—model validation and performance evaluation (growth stage) Using the 11 selected microbial markers as input features, a LightGBM-based classification model was constructed to predict the corn growth stage of the sample. The leave-one-out method was used, with one sample used as the test set each time and the remaining samples used for training. A 10-fold cross validation was performed to determine the optimal number of iterations. The model was then trained at this number of iterations to predict the remaining samples. The prediction results for each sample were saved, including whether the prediction was correct, the actual label, and the predicted label.

[0052] Figure 8 This is a confusion matrix diagram for predicting phenotypes of corn growth stages based on LightGBM. The confusion matrix shows the corresponding relationship between the prediction results and actual values ​​of the classification model in the four periods. The rows of the matrix represent the predicted values ​​of the model, the columns represent the actual values, and the numbers on the diagonal are the number of samples predicted correctly. Overall, the classification accuracy of the model in the four periods is high, but the error is mainly concentrated in the classification confusion between adjacent periods. This shows that although the model can effectively distinguish different periods, the feature similarity of adjacent periods may have a certain impact on the classification performance.

[0053] To comprehensively evaluate the classification performance of the LightGBM model, Fig. 9 As shown in the figure, we plotted the ROC curve and listed the main classification indicators. By adjusting the prediction threshold, different true positive rates TPR and false positive rates FPR were obtained, and the ROC curve was generated. The AUC value is close to 1, indicating that the model has a high discriminatory power in distinguishing different growth stages. The large AUC area indicates that the LightGBM model performs well in distinguishing differences in corn growth stages and has a good prediction effect.

[0054] Table 8 Confusion matrix parameter table of growth phase validation model Table 8 is a confusion matrix parameter table of the validation model in the growth stage (recall rate indicates the proportion of correctly identifying the category. Specificity indicates the proportion of correctly excluding non-category. Positive prediction value indicates the proportion of samples predicted as the category that are actually correct. Negative prediction value indicates the proportion of samples predicted as non-category that are actually correct. Balanced accuracy indicates the average of recall rate and specificity, which measures the overall prediction effect.), as can be seen from Table 8, the characteristics between the fourth week (DAS28) and the sixth week (DAS42) are relatively similar, which leads to some confusion in the classification of the model in these stages, and the accuracy and positive prediction value are relatively low. The eighth week (DAS56) showed relatively unique characteristics, with a high model recognition rate, but overlapped with the tenth week in some cases. The characteristics of the tenth week (DAS70) are quite different from those of other weeks, and the model has the best recognition effect in this stage, with high sensitivity and positive prediction value. These results show that as the growth stage progresses, the characteristic differences of the microbial community become increasingly obvious, especially in the tenth week, which is conducive to model differentiation.

[0055] To more comprehensively evaluate the model performance, we calculated common classification metrics and compared them with models under different feature combinations, including the following metrics: Accuracy: The overall prediction accuracy of the model for all samples, that is, the proportion of correctly classified samples to the total samples. Precision: The proportion of samples correctly predicted by the model as a certain category to all samples predicted as that category. The result is the average precision of each category. Recall: The proportion of samples correctly predicted by the model as a certain category to all samples that actually belong to that category. The result is the average recall of each category. F1 Score: The harmonic mean of precision and recall, used to comprehensively measure the performance of the model. The result is the average F1 score of each category. Specificity: The proportion of samples correctly predicted by the model as negative classes, that is, the proportion of samples in non-target classes that are correctly classified. The result is the average specificity of each category. Kappa coefficient: Measures the consistency between the model classification results and the random classification results, removing the influence of random factors, ranging from -1 to 1, and the higher the value, the better the consistency.

[0056] For multi-classification problems, the indicators (such as precision and recall) of each category are different. In this example, each category is equally important, and taking the average value can be simplified into an overall indicator, which reflects the overall performance of each indicator on all categories.

[0057] The average precision and average recall are close, indicating that the overall classification ability of the model is balanced in multi-classification tasks.

[0058] Table 9 Model evaluation table (growth stage) As shown in Table 9, "all genera" contains 2025 genera obtained by annotation; "random genera" is 75 randomly selected genera; "LDA>2" is 75 differential microbial genera with LDA scores greater than 2 screened out by LEfSe; "LDA>2.5" is 38 genera with LDA scores greater than 2.5 screened out by LEfSe; "LDA>3" is 17 genera with LDA scores greater than 3 screened out by LEfSe; the final "microbial markers" are 11 key genera obtained based on LightGBM feature screening. The model of the microbial marker feature set performed well in all indicators and was significantly better than other feature selection methods. This shows that the key microbial features screened based on LightGBM can effectively predict different growth stages and have high classification performance. In contrast, the model containing all genera and the model of randomly selected genera have lower accuracy and Kappa coefficients (0.671 and 0.561; 0.529 and 0.371, respectively), indicating that it is difficult for unscreened feature sets to accurately distinguish growth stages. Among other feature selection methods, the feature set with an LDA score greater than 2.5 performed relatively well in terms of precision (0.704) and recall (0.702), but was still lower than the performance of the "microbial marker" feature set (precision was 0.789 and recall was 0.785). This shows that the use of important features after screening can not only improve the classification performance of the model, but also significantly reduce the interference of redundant features on the prediction results. In addition, we found that certain growth stages (such as the late stage) are more easily accurately predicted by the model, which may reflect that the microbial community structure in these stages is more stable and has a higher consistency, thereby improving the accuracy and reliability of classification.

[0059] Finally, it should be noted that the above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art can still modify the technical solutions described in the aforementioned embodiments or replace some of the technical features therein by equivalents. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present invention shall be included in the protection scope of the present invention.

Claims

1. A method for screening corn microbial markers and predicting phenotypes, characterized in that: The following steps are involved: S1. Microbial sample pretreatment: Obtain corn microbial sample data, perform data cleaning and standardization operations, generate a standardized microbial abundance table, perform Alpha diversity analysis to understand the internal structure of the microbial community, and verify whether the pretreatment method is reasonable; S2.LEfSe analysis to screen differential microorganisms: group the samples in the microbial abundance table; use the Kruskal-Wallis rank sum test to screen the microbial features with significant differences; use the Wilcoxon rank sum test to compare the two groups to screen the microbial features with stable and significant differences between the groups; use linear discriminant analysis to calculate the effect size of microbial features with stable and significant differences, and finally obtain the set of microbial features with significant differences; S3. Feature importance evaluation based on LightGBM: Divide the significantly different microbial feature set screened by LEfSe into a training set and a test set; construct a LightGBM model; the LightGBM model performs computational evaluation on the importance of microbial features; screen important features to obtain microbial markers for in-depth analysis of key microbial community characteristics of maize phenotypic changes; S4. Model validation and performance evaluation: The predictive validity of the screened microbial features was verified through model training, model validation and model performance evaluation.

2. A method for screening corn microbial markers and predicting phenotypes according to claim 1, characterized in that: Step S1 analyzes the microbial community in the corn microbial sample data by high-throughput sequencing technology to obtain a microbial abundance table, wherein the high-throughput sequencing technology includes 16S rRNA gene sequencing and metagenomic sequencing.

3. A method for screening corn microbial markers and predicting phenotypes according to claim 2, characterized in that: Step S2 includes the following steps: S21. Input data and grouping: Group the corn samples in the microbial abundance table according to the experimental design. Each sample has a corresponding classification label to distinguish the differences between the groups and ensure that the subsequent analysis can accurately capture the differences in microbial characteristics of different groups; the classification label is the corn phenotype, which includes genotype and growth stage; S22. Kruskal-Wallis rank sum test to screen significant microbial features: LEfSe first performed a Kruskal-Wallis rank sum test on the microbial abundance table; The significance level was set to p value less than 0.05, and the microbial characteristics with statistically significant differences between different groups were screened out. The formula is as follows: ， in, H Represents the test statistic, which is used to determine whether there is a significant difference between groups. N represents the total number of samples, m Indicates the number of groups, n i represents the number of samples in the i-th group, R i represents the rank sum of the i-th group; S23. Wilcoxon rank sum test was used to compare the two groups to screen out the microbial characteristics with stable and significant differences between the groups: For the microbial characteristics with significant differences through Kruskal-Wallis test, Wilcoxon rank sum test was further performed to evaluate the specific differences between each group; S24. Linear discriminant analysis was used to calculate the effect size of microbial features with stable and significant differences: Linear discriminant analysis (LDA) was used to calculate the effect size of microbial features with stable and significant differences. By setting the effect size threshold, the microbial features that contributed more to the differences between groups were screened out, and then a set of microbial features with significant differences was obtained. The effect size calculation formula of LDA is: ， Among them, LDA Score represents the effect size of microbial characteristics. The larger the value, the greater the contribution of this microorganism to the difference. μ g1 , μ g2 They represent the characteristic means of the first group and the second group respectively, and σ represents the overall standard deviation of the characteristic.

4. A method for screening corn microbial markers and predicting phenotypes according to claim 3, characterized in that: Step S3 includes the following steps: S31. Data preparation: The significantly different microbial feature set screened by LEfSe is used as the input feature of the LightGBM model, and the corn phenotype is used as the prediction target of the model; the significantly different microbial feature set screened by LEfSe is divided into a training set and a test set; S32. Construct LightGBM model: Use the histogram-based decision tree construction method to discretize the abundance values ​​of microorganisms in the microbial abundance table into several intervals and construct a histogram. The model splits by selecting the feature that produces the largest reduction in the loss function, generates binary decision trees layer by layer, and finally obtains the LightGBM model. Set the basic parameters of the LightGBM model, including the number of leaf nodes num_leaves, learning rate learning_rate, and maximum depth max_depth of the tree for each tree; S33. Feature importance evaluation: The LightGBM model provides a variety of feature importance evaluation methods, including feature importance calculation based on gain, coverage, and split times. The significant difference microbial feature set is evaluated through feature importance evaluation to obtain the importance score of the significant difference microbial feature. S34. Screening important features: sorting the significantly different microbial features according to the importance scores of the significantly different microbial features obtained in step S33, and screening important features as microbial markers; S35. Predict corn phenotype based on the microbial markers obtained in step S34 to obtain the final predicted phenotype.

5. The method for screening corn microbial markers and predicting phenotypes according to claim 1, characterized in that: Step S4 includes the following steps: S41: Model training: Use the significantly different microbial features in the training set as the input features of the LightGBM model, take the corn phenotype as the prediction target of the model, train the LightGBM model on the training set, learn the association pattern between the significantly different microbial features and phenotypic changes, and finally obtain a trained LightGBM model; S42: Model validation: Validate the trained LightGBM model on the test set to evaluate the prediction performance of the model on unseen data; obtain the overall performance indicators of the model through the prediction results of the test set; use the cross-validation method to further verify the robustness of the model; S43: Model performance evaluation: Combine the test set and cross-validation results to comprehensively evaluate the performance indicators of the model.

Images