Soil clayey layer thickness digital mapping method based on Boruta feature screening
By combining Boruta feature selection and QRF model, the problem of insufficient prediction accuracy of soil clay layer thickness was solved, achieving efficient and reliable prediction of soil clay layer thickness, reducing computational costs and improving model accuracy.
Patent Information
- Application Number
- CN202511421980.7
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-09-30
- Publication Date
- 2026-03-17
AI Technical Summary
Existing methods for predicting soil clay layer thickness in the Northeast Plain of China suffer from insufficient prediction accuracy due to the scarcity of soil samples and the simple terrain. It is also difficult to select appropriate environmental covariates to improve prediction accuracy. Traditional methods are time-consuming, labor-intensive, and costly.
A Boruta-based feature selection method was adopted, which combined multi-dimensional environmental covariates such as topography, climate, biology, and soil. The resolution was unified by bilinear interpolation, and the model was trained using a quantile regression forest (QRF) model. Key features were selected by combining 10-fold cross-validation and the Boruta algorithm to construct a soil clay layer thickness prediction model.
It improves the accuracy and reliability of soil clay layer thickness prediction, reduces computational load and processing costs, provides predicted values and confidence intervals to guide subsequent investigations, and enhances the efficiency of results application and transformation.
Smart Images
Figure CN121685693A_ABST
Abstract
Description
TECHNICAL FIELD
[0001] The present application relates to the technical field of geographic information, in particular to a soil claypan thickness digital mapping method based on Boruta feature screening. BACKGROUND
[0002] Soil claypan is one of the diagnostic horizons in the United States soil classification, which is located below the soil surface layer and has a significantly higher clay content than the overlying soil layer. Soils with claypan are mainly classified as luvisols in the Chinese and American soil system classification, and belong to the grey luvisol primary unit in the United Nations world soil. Claypan is commonly found in arid and semiarid regions, and is easily formed under saline-alkali soil and specific climatic conditions, which is of great significance in agricultural, horticultural and soil science research. Its thickness is an important property of luvisol, and has a strong control effect on the surface and underground soil processes.
[0003] The formation of claypan is caused by the migration and deposition of clay in soil, which can be divided into two categories: residual claypan and depositional claypan. Residual claypan is formed by the weathering of primary minerals, while depositional claypan is caused by the migration and deposition of clay. Climate and topography are the main factors affecting the formation of claypan, among which climate plays a decisive role. However, there are still many uncertainties about the formation mechanism of claypan, such as the differences in morphological characteristics of claypan under different claypan processes, and how to accurately determine the dominant type of claypan, which are still difficult problems that need to be further explored in the academic field.
[0004] The diagnostic features of claypan include: the ratio of clay content in the accumulation layer to the leaching layer is usually greater than 1.2, and the thickness of claypan is not less than 1 / 10 of the total thickness of the leaching layer and the accumulation layer. In addition, the clay film inside or near the surface of claypan usually shows directional arrangement and other characteristics. Although these features are relatively clear, the specific impact threshold on crop growth still needs further verification, especially under different climate and soil conditions, the performance and impact of claypan may be different.
[0005] There are different views on the role of claypan in the academic field. Some scholars believe that claypan can effectively improve the water retention and nutrient content of soil, and help prevent the penetration of pollutants, thereby providing better growing conditions for plants. On the contrary, some studies have pointed out that the heavy clay texture of claypan can affect the tillage performance of soil, reduce water permeability, especially in the rainy season, which can easily lead to waterlogging and affect the growth of plant roots, and even cause root rot and root rot. Therefore, how to balance the advantages and disadvantages of claypan is still an important topic in soil science research.
[0006] In the soil mapping technology, the traditional method relies on a large number of field investigations and manual drawing, which is not only time-consuming and laborious, but also high in cost. With the rapid development of digital soil mapping technology, the spatial prediction technology based on the'soil-landscape model' gradually becomes the mainstream. This method can improve the prediction accuracy of soil thickness distribution by combining topography, climate, parent material and other factors, which is superior to the traditional geostatistical method. However, due to the high spatial heterogeneity of soil thickness and the joint influence of multiple environmental factors, especially in the northeast plain area of China, due to the scarcity of soil samples and the simple terrain, the existing prediction accuracy still faces great challenges. Therefore, how to select appropriate environmental covariates to improve the prediction accuracy is still an important problem to be solved in current technology.
[0007] Therefore, the present application aims to provide a soil claypan thickness digital mapping method based on Boruta feature screening to solve the above problems. SUMMARY
[0008] The purpose of the present application is to solve the above problems, and to provide a soil claypan thickness digital mapping method based on Boruta feature screening, which uses sparse and limited field profile data to establish a reliable soil claypan thickness prediction model, which has important practical significance for subsequent soil barrier factor reduction and land management policy improvement.
[0009] In order to achieve the above purpose, the technical scheme of the present application is as follows:
[0010] The present application provides a soil claypan thickness digital mapping method based on Boruta feature screening, which comprises the following steps:
[0011] S1, based on the position of each sample point arranged in the target area, obtaining the soil profile sample data (i.e. thickness value) of the region, selecting soil profiles with claypan characteristics, integrating the thickness information of each claypan, and forming a continuous soil claypan thickness data set;
[0012] S2, obtaining each environmental covariate grid data covering the target area and related to soil claypan forming factors, and using resampling technology (bilinear interpolation method) to unify the spatial resolution of all grid data to the same scale, forming a homogeneous environmental covariate grid;
[0013] S3, based on the position of each sample point arranged in the target area, extracting the value of each environmental covariate to the sample point, constructing the combination of sample data and each soil-forming environmental covariate factor, and forming the basis of environmental covariate feature screening;
[0014] S4, based on the overall correlation between each environmental variable and the soil claypan thickness value of each position sample point, the preliminary environmental variable subset is obtained by preliminarily removing the redundant features of the environmental covariates;
[0015] S5, based on the preliminary environmental variable subset obtained in S4, the environmental covariates corresponding to each sample point position in the subset are input into the Boruta algorithm program, and the optimal environmental variables related to the soil claypan thickness value in the environmental variable subset are determined according to the feature importance (Z-score) calculated by the Boruta algorithm program, and the environmental variable subset obtained in S4 is updated to form an optimal environmental variable group;
[0016] S6, based on the sample point positions and the profile sample thickness values of the target region, the data of each optimal environmental variable corresponding to each sample point position in the optimal environmental variable group are input into the preset machine learning model (QRF model) for model training;
[0017] In this input, the data of each optimal environmental variable corresponding to the sample point position in the optimal environmental variable group is input, and the soil claypan thickness value corresponding to the sample point position is output. The trained model (QRF model) is trained to obtain a continuous soil claypan thickness prediction model. The model verification uses 10-fold cross-validation method, and the determination coefficient (R 2 ), root mean square error (RMSE), mean absolute error (MAE) and prediction interval coverage probability (PICP) of the soil claypan thickness prediction model are obtained. Among them, the model training is repeated 50 times independently, forming 50 independent soil claypan thickness prediction models, and the average performance of the model is used as the final evaluation;
[0018] S7, based on the 50 soil claypan thickness prediction models preset in S6, the data of each optimal environmental variable corresponding to the grid sample position in the optimal environmental variable group is input, and the soil claypan thickness value corresponding to the grid sample position is output. The 5% quantile, mean, 50% quantile and 95% quantile prediction results of the soil claypan thickness value of each grid unit in the target region are gradually output. Among them, the 50 soil claypan thickness prediction models are independently repeated for prediction, and the average value of the prediction results of the 50 independent soil claypan thickness prediction models is obtained to obtain the soil claypan thickness spatial prediction map and uncertainty map covering the target region.
[0019] The specific steps of Boruta screening in step S5 are as follows:
[0020] S5A1, based on the target area layout of each sample point position and its profile sample thickness value, the preliminary environmental variable subset obtained in S4, according to the data information of each environmental variable in the preliminary environmental variable subset corresponding to the target area, pre-constructing the RF model for predicting the thickness of the soil cementation layer, and based on the pre-constructed RF model for predicting the thickness of the soil cementation layer, calculating the feature importance of each input environmental variable one by one;
[0021] S5A2, based on the target area layout of each sample point position and its profile sample thickness value, the preliminary environmental variable subset obtained in S4, generating "shadow features" consistent with the number of environmental variables in the subset, forming a feature data set to be compared;
[0022] S5A3, based on the feature importance of each input environmental variable calculated one by one by the RF prediction model for predicting the thickness of the soil cementation layer pre-constructed in S5A1, and the feature importance of the "shadow features" generated in S5A2 consistent with the number of environmental variables in the preliminary environmental variable subset obtained in S4, the feature importance is compared with each other, and then the important features and unimportant features are divided;
[0023] S5A4, based on the steps of S5A1-S5A3, under the RF model framework of Boruta algorithm, a plurality of iterations are carried out, and each input environmental variable in the RF prediction model for predicting the thickness of the soil cementation layer pre-constructed in S5A1 and the "shadow features" generated in S5A2 consistent with the number of environmental variables in the preliminary environmental variable subset obtained in S4 are combined and a new RF prediction model for predicting the thickness of the soil cementation layer is pre-constructed; based on the comparison method of feature importance in S5A3, the important features and unimportant features are re-divided under the newly pre-constructed RF prediction model for predicting the thickness of the soil cementation layer, until the important features and unimportant features divided by the step do not change significantly, then all the important features are determined as the optimal environmental variable group.
[0024] The specific steps of model training (QRF model) in step S6 are as follows:
[0025] S6A1, based on the target area layout of each sample point position and its profile sample thickness value, the data of each sample point position corresponding to each optimal environmental variable in the optimal environmental variable group is input into the pre-set QRF model for model training; in this input, the data of the sample point position corresponding to each optimal environmental variable in the optimal environmental variable group is input, and the thickness value of the soil cementation layer corresponding to the sample point position is output;
[0026] S6A2, under the QRF model framework, the model parameter optimization adopts the grid search method, that is, each parameter is searched in the parameter optimization grid to form a hyperparameter combination, the model performance under each hyperparameter combination is evaluated, the optimal hyperparameter combination is selected with the lowest RMSE as the evaluation standard; the three important parameters of the QRF model framework are the number of features considered by the decision tree at the split node (mtry), the minimum number of samples each node must contain before splitting (nodesize) and the number of decision trees (ntree); ntree is set to the default value 500; mtry is from 1 to the maximum number of features with a step of 1; nodesize is from 3 to 10 with a step of 1;
[0027] S6A3, under the optimal hyperparameter combination, a spatial prediction model is established; wherein the model training is repeated 50 times independently, forming 50 independent soil cohesive layer thickness prediction models, and the average performance of the model is used as the final evaluation, that is, the model training part is completed.
[0028] The specific steps of 10-fold cross-validation in step S6 are as follows:
[0029] S6B1, based on the positions of each sample point and the thickness values of the profile sample points arranged in the target area, all sample points are randomly and uniformly divided into 10 subsets to form the data basis for 10-fold cross-validation;
[0030] S6B2, for the 10 subsets generated in S6B1, randomly select one subset as the validation set of the model, and the remaining 9 subsets are used as the training set of the model, and perform 10 iterations to ensure that each subset is used as a validation set once, and perform comprehensive internal cross-evaluation of the model; the precision evaluation of the final model is the mean value of the 10 iterations;
[0031] S6B3, based on the steps of S6B1-S6B2, repeat the 10-fold cross-validation of the model 50 times, take the mean value of the precision evaluation index as the final result of the model evaluation, and calculate the R 2 , RMSE, MAE and PICP evaluation indexes.
[0032] The present application breaks through the traditional single source limitation in data acquisition and integration, and constructs a spatial database of 311 effective observations by 2023-2024 field investigation in the three provinces of northeast China and integration of 98 historical profiles containing clayey horizons in the three provinces of northeast China in the 2010-2018 China Soil Series Record; meanwhile, according to the soil-landscape relationship theory and SCORPAN paradigm, 81 environmental covariates are selected from four dimensions of terrain, climate, biology and soil, including 14 terrain factors derived from Google Earth Engine 30m resolution DEM, 28 climate factors from High Resolution Mountain Environment Mapping Plan and Resource Environment Science Data Platform, 17 biological factors from Google Earth Engine Sentinel-2 image and National Qinghai-Tibet Plateau Scientific Data Center, and 22 soil factors from China High Resolution National Soil Information Grid Basic Attribute Dataset, etc., and the resolution is unified to 90m through bilinear interpolation. In the feature screening link, through the combination of Pearson correlation analysis and Boruta algorithm, 20 key features are finally determined, and the contribution of terrain and biological factors is relatively low. In the modeling and evaluation process, the quantile regression forest (QRF) model is applied to the prediction of clayey horizon thickness, the model is constructed by optimizing the decision tree splitting feature number (mtry) and the minimum sample size of the node (nodesize), the “50 times repeated 10-fold cross-validation” method is adopted, and the model accuracy is evaluated by combining the coefficient of determination (R²), root mean square error (RMSE) and mean absolute error (MAE), and the reliability of uncertainty is evaluated by introducing the prediction interval coverage rate (PICP). In the uncertainty and interpretability analysis, the uncertainty index is calculated based on the QRF output quantile and the uncertainty distribution map is drawn, and the nonlinear relationship between the covariates and the clayey horizon thickness is revealed by means of feature importance sorting and accumulated local effect (ALE) map, which breaks through the limitation of the “black box” model.
[0033] Compared with the prior art, the present application has the following advantages:
[0034] The QRF model of the present application can provide a predictive value condition distribution and quantify uncertainty, the PICP of 87.1% is close to the target value of 90%, the user can be provided with “thickness result + confidence interval”, the uncertainty map can also guide subsequent supplementary investigation, reduce blind sampling cost, and improve the application transformation efficiency of the results; the feature screening system eliminates redundant variables, reduces the calculation amount and improves the operation efficiency by simplifying 81 covariates to 20, and indirectly improves the model accuracy by focusing on the core driving factors; the method clearly determines the variable priority in the large-scale plain area with gentle terrain, and can be reused for the prediction of soil layer thickness or unconventional soil properties in other large areas and sparse samples by combining double feature screening, QRF modeling and multidimensional evaluation, thereby providing a reference for similar research. BRIEF DESCRIPTION OF DRAWINGS
[0035] Figure 1 is a spatial distribution diagram of the study area, soil type and profile sample points in the embodiment of the present application;
[0036] Figure 2 is a frequency histogram of the thickness of the soil claypan layer obtained by field investigation in the embodiment of the present application;
[0037] Figure 3 is a Pearson correlation coefficient diagram between the thickness of the soil claypan layer and the environmental covariates (P<0.01) in the embodiment of the present application;
[0038] Figure 4 is an environmental covariate optimization screening diagram of the thickness of the soil claypan layer based on the Boruta algorithm in the embodiment of the present application;
[0039] Figure 5 is a model performance evaluation diagram of 50 times of 10-fold cross-validation, including R2, RMSE, MAE and PICP in the embodiment of the present application;
[0040] Figure 6 is a diagram of the predicted thickness of the soil claypan layer by the QRF model of Heilongjiang, Jilin and Liaoning in the embodiment of the present application, including the mean, median, low estimation limit and high estimation limit;
[0041] Figure 7 is a diagram of the uncertainty map of the predicted thickness of the soil claypan layer generated by the QRF model of Heilongjiang, Jilin and Liaoning in the embodiment of the present application;
[0042] Figure 8 is a diagram of the importance ranking of the environmental covariates generated by the QRF model in the embodiment of the present application;
[0043] Figure 9 is a cumulative local effect diagram of the first nine important covariates generated by the QRF model in the embodiment of the present application, and the Y axis is the thickness of the claypan layer. DETAILED DESCRIPTION
[0044] In order for those skilled in the art to better understand the present application, the technical solutions of the present application will be further described in detail below in combination with the embodiments of the present application and the accompanying drawings. Obviously, the described embodiments are only a part of the embodiments of the present application, not all. Based on the embodiments in the present application, all other embodiments obtained by those skilled in the art without creative labor should be within the scope of protection of the present application.
[0045] It should be noted that the embodiments in the present application and the features in the embodiments can be combined with each other without conflict. The present application will be described in detail below in combination with the embodiments.
[0046] Example: Taking the three northeastern provinces of Heilongjiang, Jilin, and Liaoning as examples.
[0047] 1. Materials and Methods
[0048] 1.1 Study Area
[0049] The study area was selected from the three provinces of Heilongjiang, Jilin, and Liaoning in Northeast China. This region is mostly located in the temperate zone, with smaller portions in the cold temperate and warm temperate zones. Humid and semi-humid areas coexist, and the terrain is predominantly plains and mountains. The region has a rich variety of soil types, including black soil, chernozem, meadow soil, and brown soil. The clayey layer is mainly found in leached soils. The distribution of other soil types was confirmed based on soil survey data combined with geographic coordinates. Detailed distribution is as follows... Figure 1 As shown.
[0050] 1.2 Soil Data
[0051] In 2023 and 2024, soil surveys were conducted based on data from the Second National Soil Census and auxiliary information such as elevation and NDVI. Survey areas were selected based on road accessibility, with survey points primarily set up in low-lying and high-lying areas. Soil profiles 1.5-2 meters deep were excavated at each survey point, and coordinates were recorded using a handheld GPS device. The determination of the clay layer strictly adhered to six criteria, with the thickness of the clay layer visually read from the center of the profile using a measuring tape. A total of 311 soil sample points were recorded.
[0052] 1.3 Environmental covariates
[0053] According to the soil-landscape relationship theory, soil formation follows this equation:
[0054]
[0055] In the formula, Soil represents a certain soil property; (Soils) refers to other soil properties; (Climate) refers to climatic factors; (Organism) refers to biological factors; (Relief) refers to topographic factors; (Parent material) refers to the parent material factor; (Age) is a time factor; (Geographic position) refers to spatial location factors; This study examines the interactions among these factors. Environmental covariates are selected from four aspects: topography, climate, organisms, and soil factors. These environmental factors fall into four main categories: topography, climate, organisms, and soil.
[0056] 1.3.1 Topography
[0057] The DEM data for the terrain factors ('NASA / NASADEM_HGT / 001') originated from Google Earth Engine (GEE) with an initial resolution of 30m. The terrain-derived factors were calculated using ArcMap 10.8 and include 14 terrain factors such as elevation.
[0058] 1.3.2 Climate
[0059] The climate factor variables come from two sources. The first is data products from the High-Resolution Mountain Environment Mapping Project of the Institute of Mountain Hazards and Environment, Chinese Academy of Sciences, including 23 climate factors such as MAP, with a mean of 1990s-2010s and an initial resolution of 30m; the second is from the Resource and Environmental Science Data Platform of the Chinese Academy of Sciences, including 5 climate factors such as MAE, with a mean of 1960s-2010s and an initial resolution of 1000m.
[0060] 1.3.3 Organisms
[0061] Sentinel 2 bands were obtained from the GEE platform as biological covariates, with the time range consistent with the sampling period (2023-2025), and median composites of the image bands were performed. In addition, spectral index data (mean values from 2000s to 2010s) were collected from the National Tibetan Plateau Scientific Data Center, including five biological factors such as composite NDVI maximum values, and were used together with the Sentinel 2 image bands for feature screening.
[0062] 1.3.4 Soil
[0063] Soil factor variables use soil surface property data. In Northeast China, leaching and deposition processes produce a relatively stable clay layer. Although the depth of clay migration cannot be determined, soil surface properties are likely to change as a result. This process usually takes place under acidic conditions, and the reduction of clay particles may lead to sandification of the surface soil (initial resolution 90m).
[0064] 1.3.5 Unification of Spatial Resolution of Environmental Covariates
[0065] All environmental covariates were uniformly resampled to a 90m resolution before spatial prediction modeling and mapping were performed.
[0066] 1.4 Feature Selection
[0067] 1.4.1 Pearson Correlation Analysis
[0068] Pearson correlation analysis is used to assess the degree of linear correlation between multiple pairs of continuous variables, with values ranging from -1 to 1: values close to -1 indicate a strong negative correlation, values close to 0 indicate a weak correlation, and values close to 1 indicate a strong positive correlation. This application uses IBM Statistics SPSS 19.0 software to analyze the correlation of environmental covariates, eliminating variables with insignificant correlation (P < 0.01).
[0069] 1.4.2 Boruta
[0070] Boruta is a feature selection method based on random forests. It determines a representative feature subset by constructing "shadow features" and comparing their importance with the real features. The steps are: (1) construct a random forest model and calculate feature importance; (2) create "shadow features" for comparison; (3) determine important and unimportant features based on the comparison results; (4) iterate through multiple rounds, using the important features from the previous round and "shadow features" to construct a new model, and update the set of important features until there is no significant change. It is implemented using the Boruta package in the R language, with important features as the modeling input.
[0071] 1.5 Spatial Prediction Modeling and Model Performance Evaluation
[0072] A quantile regression forest model was constructed to spatially predict the thickness of the soil clay layer, using an algorithm developed by Meinshausen. Model calibration considered three parameters: mtry, nodesize, and ntree. ntree was set to the default value of 500. The mesh settings were optimized, with mtry ranging from 1 to the maximum number of features and a step size of 1, and nodesize ranging from 3 to 10 and a step size of 1. Modeling was performed using the optimal parameters.
[0073] The model performance was evaluated using 10-fold cross-validation. The sample points were randomly divided into 10 subsets: one subset was used as the validation set, and the other nine subsets were used as the training set. This process was repeated 10 times, and the average value was used to evaluate accuracy. The model was then subjected to 50 cycles of 10-fold cross-validation, and the average accuracy metric was used as the final result. Four evaluation metrics were selected: R², RMSE, MAE, and PICP.
[0074] 1.6 Uncertainty Estimation and Model Interpretability
[0075] The soil clay layer thickness at the 0.05, 0.50, and 0.95 quantiles was predicted using a quantile regression forest model, and the uncertainty index was calculated by comparing the difference between the 0.05 and 0.95 quantiles with the median.
[0076] The dplyr package in R was used to rank the environmental covariates of the QRF model by feature importance to determine the most influential environmental variables. The ALEPlot package in R was used to generate a cumulative local effects (ALE) plot to visualize the influence of the predictor variables on the prediction of the adhesive layer thickness.
[0077] 2. Conclusion
[0078] 2.1 Descriptive Statistical Analysis
[0079] Figure 2 The thickness data of soil samples containing the clay layer are presented, ranging from 7 to 138 cm, with an average of 53.77 cm, a median of 47 cm, and a standard deviation of 28.92 cm. The skewness coefficient is 0.70 (>0.5), the kurtosis coefficient is -0.18 (<0), and the coefficient of variation is 0.54. Various mathematical transformations failed to significantly improve the model's predictive performance, and the original data was ultimately used for modeling.
[0080] 2.2 Feature Filtering Results
[0081] 2.2.1 Results of Pearson Correlation Analysis
[0082] from Figure 3 It was found that the thickness of the soil clay layer was negatively correlated with SOCD (P < 0.01, correlation coefficient -0.17 to -0.39) and positively correlated with ST (P < 0.01, correlation coefficient 0.16 to 0.45), with low correlation coefficients (|r| < 0.50). The Boruta algorithm was used to optimize feature selection using these covariates.
[0083] 2.2.2 Boruta screening results
[0084] Based on the results of the Boruta algorithm ( Figure 4 All important features were used as modeling inputs for predicting soil clay layer thickness. Soil factors included 7 variables such as ST; climate factors included 13 variables such as LMTmin.
[0085] 2.3 Model Performance
[0086] Figure 5 The results of 10-fold cross-validation over 50 iterations are presented. The upper and lower limits of four model performance evaluation metrics, including R², show small fluctuations, reducing the uncertainty caused by the randomness of data partitioning and making the evaluation results more robust. The four metrics roughly follow a normal distribution, and the mean values indicate that the model accuracy is nearly unbiased. The mean values, including R², show that the model explains 32% of the spatial variation, which is better than the results reported by D. Howell et al. and also better than most studies. This demonstrates that the constructed QRF model is reliable and can characterize spatial variation features.
[0087] 2.4 Mapping of Prediction Results and Assessment of Uncertainty
[0088] Predicted results of soil clay layer thickness in Northeast China are as follows: Figure 6As shown, the soil clay layer was thicker in the western and southwestern parts of the study area, and thinner in the northern, eastern, and southeastern parts, decreasing overall from southwest to northeast. The clay layer was relatively thicker in parts of Heilongjiang, Jilin, and Liaoning, and relatively thinner in other areas, revealing the spatial differentiation characteristics of the soil clay layer thickness.
[0089] The 0.05 and 0.95 quantiles indicate prediction uncertainty, which may be caused by factors such as large prediction intervals, high data variability, and poor model performance. The mean prediction interval coverage (PICP) is 87.1%, indicating that 90% of the prediction intervals (PI) are effective and reliable. Uncertainty map ( Figure 7 The results show that uncertainty is high in mountainous and hilly areas, and low in areas such as the Northeast Plain.
[0090] 2.5 Feature Importance and Cumulative Local Effects
[0091] 2.5.1 Feature Importance
[0092] Figure 8 The importance ranking of environmental covariates used in constructing the QRF model is presented. The results show that soil thickness (ST) is the most important covariate in the modeling, followed by AP, pH, and LMTmin. The QRF model demonstrates the importance of these features (…). Figure 8 The results of the correlation analysis with Pearson () Figure 3 The results are not entirely consistent, which suggests that the relationship is likely non-linear.
[0093] 2.5.2 Cumulative Local Effects
[0094] Figure 9 The ALE plots of the top nine key predictors using the QRF model are shown, confirming that the relationship between soil clay layer thickness and environmental covariates is non-linear. Figure 9 The thicker clay layer in soil is mainly affected by high ST, high LMTmin, and high MATmin values. The clay layer is thicker in acidic soils, and high or low AP values are conducive to the development of the clay layer. The influence of other covariates fluctuates.
[0095] The above specific embodiments are merely explanations of the present invention and are not intended to limit the present invention. After reading this specification, those skilled in the art can make modifications to these embodiments without contributing any inventive step, but as long as they are within the scope of the claims of the present invention, they are protected by patent law.
Claims
1. A soil claypan thickness digital mapping method based on Boruta feature screening, characterized in that: The method comprises the following steps: S1, based on the target area, arranging each sample point position, obtaining the soil profile sample data of the target area, selecting the soil profile with claypan characteristics from the target area, and integrating the thickness information of each claypan to form a continuous soil claypan thickness data set; S2, obtaining each environmental covariate grid data covering the target area and related to the soil claypan forming factors, and using resampling technology to unify the spatial resolution of the environmental covariate grid data to the same scale to obtain a homogeneous environmental covariate grid; S3, each sample point position corresponds to each environmental variable data in the homogeneous environmental covariate grid, and a combination of sample point data and each soil-forming environmental covariate factor is constructed to form the basis of environmental covariate feature screening; S4, based on the overall correlation between each environmental variable and the thickness value of the soil claypan at each position sample point, the preliminary environmental variable subset is formed by preliminarily removing the redundant features of the environmental variables corresponding to each sample point position obtained in S3; S5, input each environmental covariate corresponding to each sample point position in the environmental variable subset into the Boruta algorithm program, determine the optimal environmental variable related to the thickness value of the soil claypan in the environmental variable subset according to the feature importance calculated by the Boruta algorithm program, update the environmental variable subset obtained in S4, and obtain the optimal environmental variable group; S6, based on the sample point position and its profile sample thickness value arranged in the target area, the data of each optimal environmental variable in the optimal environmental variable group corresponding to each sample point position is input into the preset machine learning model for model training; In this input, the data of each optimal environmental variable in the optimal environmental variable group corresponding to the sample point position is input, and the thickness value of the soil claypan corresponding to the sample point position is output, the machine learning model is trained, a continuous soil claypan thickness prediction model is obtained, the model is verified by 10-fold cross-validation method, and the determination coefficient, root mean square error, mean absolute error and prediction interval coverage probability of the soil claypan thickness prediction model are obtained; The machine learning model is repeated 50 times independently to form 50 independent soil claypan thickness prediction models, and the average performance of the soil claypan thickness prediction model is used as the final evaluation. S7, based on the 50 soil claypan thickness prediction models preset in S6, taking the data of each optimal environmental variable in the optimal environmental variable group corresponding to the grid sample position as input, and taking the soil claypan thickness value corresponding to the grid sample position as output, gradually outputting the 5% quantile, mean, 50% quantile and 95% quantile prediction results of the soil claypan thickness value of each grid unit in the target area; wherein, the 50 soil claypan thickness prediction models are independently repeated for prediction, and the average value of the prediction results of the 50 independent soil claypan thickness prediction models is taken to obtain the soil claypan thickness spatial prediction map and the uncertainty map covering the target area.
2. The method of soil claypan thickness digital mapping based on Boruta feature selection as claimed in claim 1, wherein: The specific steps of performing Boruta screening in the step S5 are: S5A1, based on the sample point positions and the soil profile sample thickness values of the target area, pre-constructing a RF model for predicting the soil claypan thickness according to the data information of each environmental variable in the preliminary environmental variable subset obtained in S4 corresponding to the target area, and calculating the feature importance of each input environmental variable based on the pre-constructed RF prediction model for the soil claypan thickness; S5A2, based on the sample point positions and the profile sample thickness values of the target area, generating "shadow features" consistent with the number of environmental variables in the preliminary environmental variable subset obtained in S4 from the preliminary environmental variable subset obtained in S4, to form a feature data set to be compared; S5A3, comparing the feature importance of each input environmental variable calculated based on the pre-constructed RF prediction model for the soil claypan thickness in S5A1 with the feature importance of the "shadow features" consistent with the number of environmental variables in the preliminary environmental variable subset obtained in S4 generated in S5A2, and then dividing the important features and the unimportant features; S5A4, based on the steps S5A1-S5A3, performing multiple iterations under the RF model framework of the Boruta algorithm, combining each input environmental variable in the RF prediction model for the soil claypan thickness pre-constructed in S5A1 and the "shadow features" consistent with the number of environmental variables in the preliminary environmental variable subset obtained in S4 generated in S5A2, and pre-constructing a new RF prediction model for the soil claypan thickness; based on the comparison method of the feature importance described in S5A3, re-dividing the important features and the unimportant features under the newly pre-constructed RF prediction model for the soil claypan thickness until there is no significant change in the important features and the unimportant features divided in the step, and then determining all the important features as the optimal environmental variable group.
3. The method of soil claypan thickness digital mapping based on Boruta feature selection as claimed in claim 1, wherein: The specific steps of model training in the step S6 are: S6A1, based on the sample point positions and the profile sample thickness values of the target area, inputting the data of each sample point position corresponding to each optimal environmental variable in the optimal environmental variable group into the preset QRF model for model training; in this input, taking the data of each optimal environmental variable in the optimal environmental variable group corresponding to the sample point position as input, and taking the soil claypan thickness value corresponding to the sample point position as output; S6A2, under the framework of QRF model, the grid search method is used to optimize the model parameters, that is, each parameter is searched in the parameter optimization grid to form a hyperparameter combination, and the model performance under each hyperparameter combination is evaluated, and the optimal hyperparameter combination is selected with the lowest RMSE as the evaluation standard; the three important parameters of the QRF model framework are the number of features considered by the decision tree at the split node, the minimum number of samples each node must contain before splitting, and the number of decision trees; the number of decision trees is set to the default value of 500; the number of features is from 1 to the maximum number of features with a step of 1; the minimum number of samples is from 3 to 10 with a step of 1; S6A3, under the optimal hyperparameter combination, a spatial prediction model is established; wherein the model training is repeated 50 times independently, forming 50 independent soil clay layer thickness prediction models, and the average performance of the model is used as the final evaluation, that is, the model training part is completed.
4. The method of soil claypan thickness digital mapping based on Boruta feature selection as claimed in claim 3, wherein: The specific steps of the step S6 of 10-fold cross-validation are: S6B1, based on the positions of the sample points arranged in the target area and the thickness values of the profile sample points, all the sample points are randomly and uniformly divided into 10 subsets to form the data basis for 10-fold cross-validation; S6B2, for the 10 subsets generated in S6B1, randomly select one subset as the validation set of the model, and the remaining 9 subsets are used as the training set of the model, and perform 10 rounds of iteration to ensure that each subset is used as a validation set once, and perform comprehensive internal cross-evaluation of the model; the final model precision evaluation is the mean value of 10 rounds of iteration; S6B3, based on the steps of S6B1-S6B2, repeat the 10-fold cross-validation of the model 50 times, take the mean value of the precision evaluation index as the final result of the model evaluation, and calculate the determination coefficient, root mean square error, mean absolute error and prediction interval coverage probability.