Wide-area multi-level soil attribute fine mapping method based on spatial downscaling

Abstract

The invention provides a wide-area multi-level soil attribute fine mapping method based on spatial downscaling. According to the method, long-time-sequence optical data information and ground elevation information are collected through a satellite sensor, long-time-sequence temperature and rainfall data are collected through a meteorological station, and spatial position auxiliary data extracted from the long-time-sequence temperature and rainfall data are combined with existing national-level soil grid database data to construct a high-quality data set; the method is used for training a random forest model constructed based on an integrated learning technology so as to realize high-precision and high-reliability refined soil attribute mapping under a wide-area scale. The method is reasonable in conception, high-precision and high-reliability refined soil attribute mapping under a wide-area scale is realized, a specific region multi-level soil attribute database with higher resolution can be constructed, and the method is suitable for popularization and application.

Description

Technical Field

[0001] The present invention relates to the technical field of soil property database construction, and in particular to a wide-area multi-level soil property fine mapping method based on spatial downscaling. Background Art

[0002] Accurate information on the spatial distribution of soil properties (such as soil particle size, organic carbon content, soil bulk density, etc.) is essential for accurately understanding environmental issues such as ecosystem services, food security, and heavy metal pollution. Due to the increasing demand for quantitative soil information at a wide scale to solve these problems, and the shortcomings of traditional soil survey methods in complex environments that are difficult to reach, time-consuming and costly, the study of digital soil mapping has become popular worldwide. In the past few decades, the main method of soil property mapping has been to formalize the empirical quantitative description of soil formation, which can be represented by constructing models of soil properties and environmental factors (including topography, vegetation, soil parent material, climate, etc.). Based on these technologies, various countries have established a series of soil digital databases. However, except for a few regions that have built a 30-meter resolution soil property database based on dense soil profile data, most of the existing soil databases are built in response to the "GlobalSoilMap" plan proposed by the International Soil Organization to draw national-level 100-meter spatial resolution soil property maps. Larger global databases, such as SoilGrids and the Harmonized World Soil Database (HWSD), use coarser resolutions. The low spatial resolution of these data limits their refined applications at the regional scale. For example, effective soil products required for precision agriculture may have a spatial resolution of at least 30 meters. Other fields that focus on soil data, such as vehicle off-road path planning and geotechnical engineering, may require higher resolution soil information. The main obstacle to building a more refined soil property database can be attributed to the severe lack of soil profile data, because data-driven soil property mapping uses a large amount of legacy soil profile data collected since the early 20th century to expand the sample size to achieve higher fitting accuracy. However, due to the poor positioning technology in the past, position errors of more than 50 meters may result, forcing various research teams to reduce mapping resolution to reduce uncertainty.

[0003] One way to obtain a finer soil map without additional samples is to spatially downscale the existing soil property map, that is, to generate a new finer soil map from the original soil map. Similarly, in order to adapt to modern raster soil digital maps, statistical downscaling methods have been introduced, which use relevant refined auxiliary variable data to parse coarse information to obtain soil property information with higher spatial resolution. At present, a foreign research team has proposed a cyclic iteration-based dissever algorithm framework based on terrain data and optical remote sensing auxiliary variables combined with a machine learning model for spatial downscaling of soil organic carbon in farmland. The algorithm uses a generalized additive model or other machine learning model to construct the relationship between fine resolution covariates and the coarse resolution original map. During the iteration process, the algorithm first uses nearest neighbor resampling to interpolate the coarse resolution map of the target variable to the grid used by the fine resolution covariate, and then regresses the fine grid value of the target variable with a series of covariates, and upsamples the predicted value of this regression model to the original base map resolution by block averaging. If the number of iterations is greater than 1, check whether the upsampled estimate has changed from the previous iteration. If the estimated change is greater than a predefined threshold, proceed to the next step, otherwise stop. Each iteration calculates the quality difference of each coarse grid pixel, i.e., the difference between the mean of the downscaled prediction and the original value of each pixel, and uses it to correct the bias factor of the fine grid estimate for guiding the next iteration.

[0004] There are two problems with the existing technology. First, the algorithm based on loop iteration consumes a lot of computing resources. The refined soil property map at a wide scale contains a large number of pixels. The large-scale repeated calculation steps in the algorithm lead to extremely low mapping efficiency. In addition, due to the complex environment and spatial patterns of soil formation at a wide scale, it is difficult to effectively characterize the spatial heterogeneity caused by regional differences in the soil formation process and deep soil characteristics by relying solely on the surface and terrain information provided by remote sensing data.

[0005] At the same time, in the existing technology, the algorithm based on loop iteration consumes a lot of computing resources, and the large-scale repeated calculation steps in the algorithm lead to extremely low mapping efficiency; due to the complex environment and spatial patterns of soil formation at a wide scale, it is difficult to effectively characterize the spatial heterogeneity caused by regional differences in the soil formation process and deep-level soil characteristics by relying solely on the surface and terrain information provided by remote sensing data.

[0006] To sum up, it is necessary to make further innovations on existing technologies. Summary of the invention

[0007] In response to the technical problems existing in the above-mentioned background technology, the present invention proposes a wide-area multi-level soil property fine-grained mapping method based on spatial downscaling. The method has a reasonable conception and realizes high-precision and high-reliability fine-grained soil property mapping at a wide-area scale, and can construct a higher-resolution multi-level soil property database for a specific region.

[0008] In order to solve the above technical problems, the present invention provides a wide-area multi-level soil property fine mapping method based on spatial downscaling, which mainly includes the following steps:

[0009] 1) Spatial downscaling based on random forest

[0010] In the spatial downscaling framework based on random forest, the spatial resolution of the soil property map to be spatially downscaled is taken as the original resolution, and the resolution of the auxiliary variables is taken as the target resolution; the image generation tool is used to generate raster images of spatial variables for the images of the original resolution and the target resolution respectively; the optical remote sensing index, elevation and terrain-related variables and temperature and rainfall variables of the target resolution are resampled to the original resolution, and the values ​​of the auxiliary variables and soil properties of each corresponding pixel are extracted in combination with the soil property maps of the spatial variables, target depth and type of the original resolution to form a corresponding relationship and construct a model training database; a random forest model is constructed for each target soil property according to the model training database, and then the number of trees, maximum depth, minimum number of samples required for internal nodes, and number of features to be considered when finding the best split of the random forest model are automatically adjusted to obtain any soil physical and chemical property map of any soil layer depth: the trained random forest model is then used to predict soil properties based on the auxiliary variables of the target resolution to obtain a soil property map of the target resolution;

[0011] 2) Residual correction based on spline interpolation

[0012] First, the trained random forest model is used to predict the target soil properties at the original resolution. Then, the difference between the original soil variables and the corresponding predicted variables is calculated to obtain the soil property residual map at the original resolution, which characterizes the unpredictable part of the target soil variables. Then, the soil residual map at the original resolution is regarded as discrete spatial points, and the spline interpolation method that is most suitable for regularly spaced data is used to spatially interpolate the soil residuals at the original resolution to obtain the soil property residual map at the target resolution. Finally, the soil property map at the target resolution predicted by the model is corrected by adding the soil property residual map at the target resolution generated by spatial interpolation to obtain the soil property map at the original spatial resolution predicted by the corrected random forest model.

[0013] The wide-area multi-level soil property fine mapping method based on spatial downscaling, wherein: in the step 1), the GDAL tool in Python is used to generate raster images of spatial variables for images of two resolutions respectively; the raster images of spatial variables include the upper left corner, lower left corner, upper right corner, lower right corner, central Euclidean distance field and longitude and latitude fields of the image.

[0014] The wide-area multi-level soil property fine-grained mapping method based on spatial downscaling, wherein: in the step 1), any one of the sampling methods of the bilinear method and the nearest neighbor method is used to resample the optical remote sensing index, elevation and terrain-related variables and temperature and rainfall variables of the target resolution to the original resolution.

[0015] The wide-area multi-level soil property fine mapping method based on spatial downscaling, wherein the specific process of constructing the model training database in step 1) is:

[0016] 1.1.1) Data preparation

[0017] Prediction variables include high-resolution long-term remote sensing images, land cover types, terrain-related variables, and long-term meteorological composite data; target variables use low-resolution national or global soil property data;

[0018] 1.1.2) Data preprocessing

[0019] Resample the target resolution predictor variables to the same low resolution as the soil data using bilinear interpolation or nearest neighbor interpolation methods depending on the type of predictor variable, so that data of different resolutions can be aligned for analysis at the same scale;

[0020] 1.1.3) Spatial variable extraction

[0021] Extract spatial variables from low-resolution pixels and merge them into predictive variables;

[0022] 1.1.4) Data Fusion

[0023] The resampled data, spatial variables and original soil attribute data are combined into a comprehensive database, with each pixel of the original resolution as a sample record; each record should contain all input variables and target output variables;

[0024] 1.1.5) Feature Selection

[0025] Pearson correlation analysis was used to screen the predictor variables to reduce the impact of multicollinearity.

[0026] The wide-area multi-level soil property fine mapping method based on spatial downscaling, wherein the specific process of automatic adjustment in step 1) is:

[0027] 1.2.1) Define the search range for hyperparameters

[0028] Define the hyperparameters to be tuned and their corresponding candidate values;

[0029] 1.2.2) Setting up GridSearchCV

[0030] Specify the number of cross-validations and define the scoring criteria to be used;

[0031] 1.2.3) Perform a grid search

[0032] Run a grid search to train the model multiple times to find the best hyperparameter combination, cross-validate all parameter combinations, and find the best performing model.

[0033] The wide-area multi-level soil property fine mapping method based on spatial downscaling, wherein the arbitrary soil physical and chemical property map of any soil layer depth in the step 1) is expressed as:

[0034]

[0035] In the above formula (1), LR represents the original low-resolution variable; Represents any soil physical and chemical property map at any soil depth; Relief LR Represents terrain-related variables; Organism & Parent LR Represents variables related to surface cover, organic life activities and soil parent material; Climate LR Represents long-term meteorological data; Spatial LR Represents a spatial variable.

[0036] The wide-area multi-level soil property fine mapping method based on spatial downscaling, wherein: the terrain-related variable Relief LR These include slope, terrain wetness index, and Melton Ruggedness Coefficient;

[0037] Variables related to surface cover, organic life activities and soil parent material Organism & Parent LR , mainly composed of long-term optical remote sensing indices, including the Normalized Difference Vegetation Index and the Clay Index;

[0038] The long-term meteorological data Climate LR , including multi-year average rainfall, multi-year average temperature, coldest quarter rainfall and wettest quarter temperature;

[0039] The spatial variable Spatial LR , including 5 Euclidean distance fields and longitude and latitude fields;

[0040] Any soil physical and chemical property map of any soil layer depth The soil physical and chemical properties can be soil organic carbon, soil bulk density, cation exchange coefficient, total phosphorus, total potassium, total nitrogen, soil gravel content, pH value, soil mechanical composition, and the soil layer depth can be from the surface to 2m underground.

[0041] The wide-area multi-level soil property fine mapping method based on spatial downscaling, wherein the soil property map of the target resolution in step 1) is expressed as:

[0042]

[0043] Wherein, HR in the above formula (2) represents the variable of target resolution; Soil property map representing the target resolution.

[0044] The wide-area multi-level soil property fine mapping method based on spatial downscaling, wherein: in the step 1), a random forest model is constructed for each target soil property separately using the RandomForestRegressor module in Python 3.9 according to the model training database, and the specific construction process is:

[0045] Assume that each tree in the random forest is an independent decision tree model. Given an input vector x = (x1, x2, ..., x p ), the decision tree divides the data into smaller subsets by binary splitting until it reaches a leaf node, which contains the prediction result h(x); the random forest makes predictions by integrating multiple decision trees, and the final prediction value is the average of the prediction values ​​of all decision trees:

[0046]

[0047] In the above formula (3), h i (x) is the prediction result of the i-th tree for sample x, and n is the total number of decision trees;

[0048] The GridSearch module is then used to automatically adjust the number of decision trees in the random forest model, the maximum depth, the minimum number of samples required for internal nodes, and the number of features to consider when finding the best split.

[0049] The wide-area multi-level soil property fine mapping method based on spatial downscaling, wherein the soil property residual map at the original resolution in step 2) is expressed as:

[0050]

[0051] Among them, in the above formula (4), LR represents the original low-resolution variable; Residual maps of soil properties representing the original resolution; represents the original soil property map that was spatially downscaled; soil property maps representing the original spatial resolution predicted by the random forest model;

[0052] The soil property residual map of the target resolution in step (2) is expressed as:

[0053]

[0054] Wherein, HR in the above formula (5) represents the variable of target resolution; Residual maps of soil properties representing the target resolution;

[0055] The original spatial resolution soil property map predicted by the corrected random forest model in step (2) is expressed as:

[0056]

[0057] Among them, in the above formula (6) Soil property maps representing the original spatial resolution predicted by the calibrated random forest model.

[0058] The wide-area multi-level soil property fine mapping method based on spatial downscaling, wherein the specific process of using the spline interpolation method most suitable for regularly spaced data to perform spatial interpolation on the soil residuals of the original resolution in step 2) is:

[0059] 2.1) Determine the interpolation area

[0060] Calculate the grid size of the high-resolution image, i.e. the number of rows and columns, according to the target resolution, and ensure that the coordinate system and spatial range of the input image are consistent during the interpolation process;

[0061] 2.2) Apply spline interpolation

[0062] Extract the pixel values ​​of the original image and their corresponding spatial coordinates and use them as the input points of the spline interpolation. Use the spline interpolation method to generate the interpolation function. Select the two-dimensional spline interpolation to adapt to the two-dimensional spatial characteristics of the raster image. The two-dimensional spline interpolation is a piecewise quadratic polynomial interpolation method, which is used to construct a series of quadratic polynomials between given data points to ensure smoothness in each interval. The goal of the two-dimensional spline interpolation is to find a set of quadratic polynomials so that the function value is continuous at each data point and the derivatives at these points are also continuous. When there are n+1 data points, the quadratic spline interpolation function S(x) is in the following form:

[0063] For each interval [xi ,x i+1 ], define a quadratic polynomial S i (x):

[0064] S i (x) = a i (xx i ) 2 +b i (xx i )+c i (7);

[0065] In the above formula (7), a i , b i 、c i is the interpolation point condition S i (x i )=y i , Continuity condition S i (x i+1 )=S i+1 (x i+1 ), derivative continuity condition S′ i (x i+1 ) = S′ i+1 (x i+1 ) and boundary conditions S″0(x0)=0 and S″ n-1 (x n )=0 determined coefficient;

[0066] 2.3) Generate high-resolution images

[0067] The interpolation results are converted into raster images and then back into raster image format. It is necessary to ensure that the new image maintains the same coordinate system and spatial reference as the target resolution image to obtain the soil property residual map of the target resolution.

[0068] By adopting the above technical solution, the present invention has the following beneficial effects:

[0069] The wide-area multi-level soil property fine-grained mapping method based on spatial downscaling in the present invention is reasonably conceived, adopts a non-iterative spatial downscaling algorithm framework to adapt to the computational efficiency requirements under wide-area scales, and further ensures the accuracy of spatial downscaling through residual correction based on spline interpolation; in addition, long-term meteorological environmental variables are introduced to provide effective information about multi-level soil formation, and spatial pseudo variables are introduced to maximize the spatial context information of the decomposed image, so that the refined image retains more spatial patterns of the original image, avoids generating unrealistic and discontinuous patterns in the data, and finally realizes the efficient wide-area multi-level soil property database fine-grained construction.

[0070] The present invention integrates long-term remote sensing data, meteorological data, and spatial data, and combines existing national soil raster data to construct a spatial downscaling framework based on random forests, thereby realizing the construction of a multi-level refined soil property database for a specific region. Therefore, the present invention has good application prospects.

[0071] The main features or advantages of the present invention are embodied in the following aspects:

[0072] (1) Wide-area multi-level soil property mapping and database construction based on spatial downscaling

[0073] Based on the existing national soil property database, a spatial downscaling framework of coarse-resolution soil properties based on ensemble learning technology was trained to achieve rapid mapping of more refined soil properties at a wide scale.

[0074] (2) Non-iterative spatial downscaling algorithm framework combined with residual correction for soil property research

[0075] In order to solve the problem that the original algorithm based on loop iteration consumes large computing resources and the repeated computing steps lead to low mapping efficiency, the present invention adopts a non-iterative spatial downscaling algorithm framework, namely, spatial downscaling based on random forest and residual correction based on spline interpolation, to meet the computing efficiency requirements under wide-area scales, and further ensures the accuracy of spatial downscaling through residual correction based on spline interpolation;

[0076] (3) Introducing long-term meteorological variables for spatial downscaling of soil properties over a wide area and at multiple levels

[0077] The synthetic results of long-term meteorological variables are introduced, that is, high spatial resolution meteorological data are added to effectively characterize the spatial heterogeneity caused by regional differences in the soil formation process and deep-level soil characteristics. It is suitable for study areas with complex soil formation environments and spatial patterns at a wide scale.

[0078] (4) Introducing spatial pseudo variables for spatial downscaling of wide-area and multi-level soil properties

[0079] Spatial pseudo variables are introduced, that is, longitude and latitude and Euclidean distance field variables are added to maximize the spatial context information of the decomposed image, so that the refined image retains more spatial patterns of the original image, avoids the generation of unrealistic and discontinuous patterns in the data, and ultimately achieves efficient and refined construction of a wide-area multi-level soil property database. BRIEF DESCRIPTION OF THE DRAWINGS

[0080] In order to more clearly illustrate the specific implementation methods of the present invention or the technical solutions in the prior art, the drawings required for use in the specific implementation methods or the description of the prior art will be briefly introduced below. Obviously, the drawings in the following description are some implementation methods of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0081] Figure 1 Flow chart of the wide-area multi-level soil property refinement mapping method based on spatial downscaling. DETAILED DESCRIPTION

[0082] The technical solution of the present invention will be described clearly and completely below in conjunction with the accompanying drawings. 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.

[0083] The present invention is further explained below in conjunction with specific implementation modes.

[0084] like Figure 1 As shown, the present embodiment provides a wide-area multi-level soil property fine-grained mapping method based on spatial downscaling, which is to construct a high-quality data set by combining the existing national soil raster database data with the long-term optical data information and ground elevation information collected by satellite sensors, the long-term temperature and rainfall data collected by meteorological stations, and the spatial position auxiliary data extracted therefrom, for training the random forest model constructed based on the ensemble learning technology, thereby realizing high-precision and high-reliability fine-grained soil property mapping at a wide-area scale, and a higher-resolution multi-level soil property data set for a specific region.

[0085] The wide-area multi-level soil property fine mapping method based on spatial downscaling of the present invention specifically includes the following steps:

[0086] 1) Spatial downscaling based on random forest

[0087] In this spatial downscaling framework, the spatial resolution of the soil property map to be spatially downscaled is taken as the original resolution, and the resolution of the auxiliary variables is taken as the target resolution.

[0088] Use the GDAL tool in Python to generate raster images of spatial variables for images of two resolutions, including the upper left corner, lower left corner, upper right corner, lower right corner, central Euclidean distance field, and longitude and latitude fields.

[0089] The optical remote sensing index, elevation and terrain-related variables, and temperature and rainfall variables of the target resolution are resampled to the original resolution through the bilinear or nearest neighbor method. Combined with the spatial variables of the original resolution and the soil property map of the target depth and type, the auxiliary variables and soil property values ​​of each corresponding pixel are extracted to form a corresponding relationship and construct the database required for model training.

[0090] According to the model training database built based on the original resolution image model, a random forest model is constructed for each target soil attribute using the RandomForestRegressor module in Python 3.9. The specific construction process is as follows:

[0091] Assume that each tree in the random forest is an independent decision tree model. Given an input vector x = (x1, x2, ..., x p ), the decision tree divides the data into smaller subsets by binary splitting until it reaches a leaf node, which contains the prediction result h(x); the random forest makes predictions by integrating multiple decision trees, and the final prediction value is the average of the prediction values ​​of all decision trees:

[0092]

[0093] In the above formula (3), h i (x) is the prediction result of the i-th tree for sample x, and n is the total number of decision trees;

[0094] Then use the GridSearch module to automatically adjust the number of decision trees, the maximum depth, the minimum number of samples required for internal nodes, and the number of features to be considered when finding the best split of the random forest model to improve the accuracy and generalization ability of the model, and obtain any soil physical and chemical property map at any soil depth:

[0095]

[0096] In the above formula (1), LR represents the original low-resolution variable. Represents any soil physical and chemical properties at any soil depth. The types can be soil organic carbon (SOC), soil bulk density (BD), cation exchange coefficient (CEC), total phosphorus (TP), total potassium (TK), total nitrogen (TN), soil gravel content (CF), pH value, soil mechanical composition (clay, silt, sand content), and the depth can be from the surface to 2m underground; Relief LR Represents terrain-related variables, including slope, terrain moisture index, Melton ruggedness coefficient, etc.; Organism & Parent LRRepresents variables related to surface cover, organic life activities and soil parent material, mainly composed of long-term optical remote sensing indices, including normalized vegetation index, clay index, etc.; Climate LR Represents long-term meteorological data, including multi-year average rainfall, multi-year average temperature, coldest quarter rainfall, wettest quarter temperature, etc.; Spatial LR Represents spatial variables, including 5 Euclidean distance fields and longitude and latitude fields.

[0097] Then use the trained random forest model to predict soil properties based on the auxiliary variables of the target resolution to obtain the soil property map values ​​of the target resolution:

[0098]

[0099] Wherein, HR in the above formula (2) represents the variable of target resolution; Represents the soil property map value at the target resolution.

[0100] The specific process of building the model training database in step 1) above is as follows:

[0101] 1.1.1) Data preparation

[0102] Prediction variables include high-resolution long-term remote sensing images of the target (including optical remote sensing indices such as NDVI and EVI), land cover type, terrain-related variables (elevation, slope, aspect, Melton ruggedness coefficient, terrain moisture index, etc.), and long-term meteorological composite data (temperature, rainfall, seasonal variation). Target variables use low-resolution national or global soil attribute data, such as soil texture composition, organic matter content, pH value, etc.

[0103] 1.1.2) Data preprocessing

[0104] Resample the target resolution predictor variables to the same low resolution as the soil data using bilinear interpolation or nearest neighbor interpolation methods depending on the type of predictor variable, so that data of different resolutions can be aligned for analysis at the same scale;

[0105] 1.1.3) Spatial variable extraction

[0106] Extract spatial variables such as latitude and longitude, Euclidean distance, etc. from low-resolution pixels and merge them into prediction variables;

[0107] 1.1.4) Data Fusion

[0108] The resampled data, spatial variables and original soil attribute data are combined into a comprehensive database, with each pixel of the original resolution as a sample record; each record should contain all input variables and the target output variable (a certain soil attribute value);

[0109] 1.1.5) Feature Selection

[0110] Pearson correlation analysis was used to screen the predictor variables to reduce the impact of multicollinearity.

[0111] The specific process of automatic adjustment in the above step (1) is as follows:

[0112] 1.2.1) Define the search range for hyperparameters

[0113] Define the hyperparameters to be adjusted and the corresponding candidate values; common hyperparameters include n_estimators (number of decision trees), max_depth (maximum depth of the tree), min_samples_split (minimum number of samples required to split a node), etc.

[0114] 1.2.2) Setting up GridSearchCV

[0115] Specify the number of cross validations (cv) and define the scoring criteria to be used (scoring);

[0116] 1.2.3) Perform a grid search

[0117] Run a grid search to train the random forest model multiple times to find the best hyperparameter combination, cross-validate all parameter combinations, and find the best performing random forest model parameter combination;

[0118] 2) Residual correction based on spline interpolation

[0119] The spatial downscaling process introduces a residual correction block to optimize the accuracy of soil properties by considering the unpredictable parts of the machine learning model. First, the trained random forest model is used to predict the target soil properties at the original resolution; then, the difference between the original soil variables and the corresponding predicted variables is calculated by the following formula (3) to obtain the soil property residual map value at the original resolution, which represents the unpredictable part of the target soil variable:

[0120]

[0121] In the above formula (4), LR represents the original low-resolution variable. Represents the soil property residual map value at the original resolution, represents the original soil property map value after spatial downscaling, Represents the original spatial resolution soil property map values ​​predicted by the random forest model.

[0122] Next, the soil residual map of the original resolution is regarded as a discrete spatial point, and the spline interpolation method that is most suitable for regularly spaced data is used to spatially interpolate the soil residual of the original resolution to obtain the soil attribute residual map value of the target resolution:

[0123]

[0124] In the above formula (5), HR represents the variable of target resolution, Represents the soil property residual map value at the target resolution.

[0125] Finally, the soil property map of the target resolution predicted by the model is compared with the soil property residual map of the target resolution generated by spatial interpolation. Add to correct:

[0126]

[0127] Among them, in the above formula (6) Represents the soil property map values ​​at the original spatial resolution predicted by the calibrated random forest model.

[0128] The specific process of spatially interpolating the soil residuals of the original resolution using the spline interpolation method that is most suitable for regularly spaced data in step 2) above is as follows:

[0129] 2.1) Determine the interpolation area

[0130] Calculate the grid size (i.e., number of rows and columns) of the high-resolution image according to the target resolution, and ensure that the coordinate system and spatial range of the input image are consistent during the interpolation process;

[0131] 2.2) Apply spline interpolation

[0132] Extract the pixel values ​​of the original image and their corresponding spatial coordinates, which will be used as the input points of the spline interpolation; use the spline interpolation method to generate the interpolation function, and select the two-dimensional spline interpolation to adapt to the two-dimensional spatial characteristics of the raster image; two-dimensional spline interpolation is a piecewise quadratic polynomial interpolation method, which is used to construct a series of quadratic polynomials between given data points to ensure smoothness in each interval; the goal of two-dimensional spline interpolation is to find a set of quadratic polynomials so that the function value is continuous at each data point and the derivatives at these points are also continuous; when there are n+1 data points, the quadratic spline interpolation function S(x) is in the following form:

[0133] For each interval [x i ,x i+1 ], define a quadratic polynomial S i (x):

[0134] Si (x) = a i (xx i ) 2 +b i (xx i )+c i (7);

[0135] In the above formula (7), a i , b i 、c i is the interpolation point condition S i (x i )=y i , Continuity condition S i (x i+1 )=S i+1 (x i+1 ), derivative continuity condition S′ i (x i+1 ) = S′ i+1 (x i+1 ) and boundary conditions S″0(x0)=0 and S″ n-1 (x n )=0 determined coefficient;

[0136] 2.3) Generate high-resolution images

[0137] The interpolation results are converted into raster images and converted into raster image format. It is necessary to ensure that the new image maintains the same coordinate system and spatial reference as the target resolution image to obtain the soil property residual map of the target resolution.

[0138] The present invention adopts a non-iterative spatial downscaling algorithm framework to adapt to the computational efficiency requirements at a wide-area scale, and further ensures the accuracy of spatial downscaling through residual correction based on spline interpolation; in addition, long-term meteorological environmental variables are introduced to provide effective information about multi-level soil formation, and spatial pseudo variables are introduced to maximize the spatial context information of the decomposed image, so that the refined image retains more spatial patterns of the original image, avoiding the generation of unrealistic and discontinuous patterns in the data, and ultimately achieving efficient and refined construction of a wide-area multi-level soil property database.

[0139] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit it. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the aforementioned embodiments, or replace some or all of the technical features therein with equivalents. However, these modifications or replacements do not cause the essence of the corresponding technical solutions to deviate from the scope of the technical solutions of the embodiments of the present invention.

Claims

1. A wide-area multi-level soil property fine-grained mapping method based on spatial downscaling, characterized in that: The main steps include: 1) Spatial downscaling based on random forest In the random forest-based spatial downscaling framework, the spatial resolution of the soil property map to be spatially downscaled is taken as the original resolution, and the resolution of the auxiliary variables is taken as the target resolution; Use the image generation tool to generate spatial variable raster images for the original resolution and target resolution images respectively; The optical remote sensing index, elevation and terrain-related variables, and temperature and rainfall variables of the target resolution are resampled to the original resolution. Combined with the spatial variables of the original resolution, the soil property map of the target depth and type, the auxiliary variables and soil property values ​​of each corresponding pixel are extracted to form a corresponding relationship and build a model training database; a random forest model is constructed for each target soil property according to the model training database, and then the number of trees, the maximum depth, the minimum number of samples required for internal nodes, and the number of features to be considered when finding the best split of the random forest model are automatically adjusted to obtain any soil physical and chemical property map of any soil layer depth: the trained random forest model is then used to predict soil properties based on the auxiliary variables of the target resolution to obtain a soil property map of the target resolution; 2) Residual correction based on spline interpolation First, the trained random forest model is used to predict the target soil properties at the original resolution. Then, the difference between the original soil variables and the corresponding predicted variables is calculated to obtain the soil property residual map at the original resolution, which characterizes the unpredictable part of the target soil variables. Then, the soil residual map at the original resolution is regarded as discrete spatial points, and the spline interpolation method that is most suitable for regularly spaced data is used to spatially interpolate the soil residuals at the original resolution to obtain the soil property residual map at the target resolution. Finally, the soil property map at the target resolution predicted by the model is corrected by adding the soil property residual map at the target resolution generated by spatial interpolation to obtain the soil property map at the original spatial resolution predicted by the corrected random forest model.

2. The wide-area multi-level soil property fine mapping method based on spatial downscaling according to claim 1, characterized in that: In the step 1), the GDAL tool in Python is used to generate raster images of spatial variables for images of two resolutions respectively; the raster images of spatial variables include the upper left corner, lower left corner, upper right corner, lower right corner, central Euclidean distance field and longitude and latitude fields of the image.

3. The wide-area multi-level soil property fine mapping method based on spatial downscaling according to claim 1, characterized in that: In the step 1), any one of the bilinear method and the nearest neighbor method is used to resample the optical remote sensing index, elevation and terrain related variables and temperature and rainfall variables of the target resolution to the original resolution.

4. The wide-area multi-level soil property fine mapping method based on spatial downscaling according to claim 1 is characterized in that: The specific process of constructing the model training database in step 1) is as follows: 1.1.1) Data preparation Prediction variables include high-resolution long-term remote sensing images, land cover types, terrain-related variables, and long-term meteorological composite data; target variables use low-resolution national or global soil property data; 1.1.2) Data preprocessing Resample the target resolution predictor variables to the same low resolution as the soil data using bilinear interpolation or nearest neighbor interpolation methods depending on the type of predictor variable, so that data of different resolutions can be aligned for analysis at the same scale; 1.1.3) Spatial variable extraction Extract spatial variables from low-resolution pixels and merge them into predictive variables; 1.1.4) Data Fusion The resampled data, spatial variables and original soil attribute data are combined into a comprehensive database, with each pixel of the original resolution as a sample record; each record should contain all input variables and target output variables; 1.1.5) Feature Selection Pearson correlation analysis was used to screen the predictor variables to reduce the impact of multicollinearity.

5. The wide-area multi-level soil property fine mapping method based on spatial downscaling according to claim 1 is characterized in that: The specific process of automatic adjustment in step 1) is as follows: 1.2.1) Define the search range for hyperparameters Define the hyperparameters to be tuned and their corresponding candidate values; 1.2.2) Setting up GridSearchCV Specify the number of cross-validations and define the scoring criteria to be used; 1.2.3) Perform a grid search Run a grid search to train the model multiple times to find the best hyperparameter combination, cross-validate all parameter combinations, and find the best performing model.

6. The wide-area multi-level soil property fine mapping method based on spatial downscaling according to claim 1 is characterized in that the arbitrary soil physical and chemical property map of any soil layer depth in step 1) is expressed as: In the above formula (1), LR represents the original low-resolution variable; Represents any soil physical and chemical property map at any soil depth; Relief LR Represents terrain-related variables; Organism & Parent LR Represents variables related to surface cover, organic life activities and soil parent material; Climate LR Represents long-term meteorological data; Spatial LR Represents a spatial variable.

7. The wide-area multi-level soil property fine mapping method based on spatial downscaling according to claim 6, characterized in that: the terrain-related variable Relief LR These include slope, terrain wetness index, and Melton Ruggedness Coefficient; Variables related to surface cover, organic life activities and soil parent material Organism & Parent LR , mainly composed of long-term optical remote sensing indices, including the Normalized Difference Vegetation Index and the Clay Index; The long-term meteorological data Climate LR , including multi-year average rainfall, multi-year average temperature, coldest quarter rainfall and wettest quarter temperature; The spatial variable Spatial LR , including 5 Euclidean distance fields and longitude and latitude fields; Any soil physical and chemical property map of any soil layer depth The soil physical and chemical properties can be soil organic carbon, soil bulk density, cation exchange coefficient, total phosphorus, total potassium, total nitrogen, soil gravel content, pH value, soil mechanical composition, and the soil layer depth can be from the surface to 2m underground.

8. The wide-area multi-level soil property fine mapping method based on spatial downscaling according to claim 1 is characterized in that: The soil property map of the target resolution in step 1) is expressed as: Wherein, HR in the above formula (2) represents the variable of target resolution; Soil property map representing the target resolution.

9. The wide-area multi-level soil property fine mapping method based on spatial downscaling according to claim 1, characterized in that: In step 1), a random forest model is constructed for each target soil property using the RandomForestRegressor module in Python 3.9 according to the model training database. The specific construction process is as follows: Assume that each tree in the random forest is an independent decision tree model. Given an input vector x = (x1, x2, ... x p ), the decision tree divides the data into smaller subsets by binary splitting until it reaches a leaf node, which contains the prediction result h(x); the random forest makes predictions by integrating multiple decision trees, and the final prediction value is the average of the prediction values ​​of all decision trees: In the above formula (3), h i (x) is the prediction result of the i-th tree for sample x, and n is the total number of decision trees; The GridSearch module is then used to automatically adjust the number of decision trees in the random forest model, the maximum depth, the minimum number of samples required for internal nodes, and the number of features to consider when finding the best split.

10. The wide-area multi-level soil property fine mapping method based on spatial downscaling according to claim 1, characterized in that: The soil property residual map at the original resolution in step 2) is expressed as: Among them, in the above formula (4), LR represents the original low-resolution variable; Residual maps of soil properties representing the original resolution; represents the original soil property map that was spatially downscaled; soil property maps representing the original spatial resolution predicted by the random forest model; The soil property residual map of the target resolution in step (2) is expressed as: Wherein, HR in the above formula (5) represents the variable of target resolution; Residual maps of soil properties representing the target resolution; The original spatial resolution soil property map predicted by the corrected random forest model in step (2) is expressed as: Among them, in the above formula (6) Soil property maps representing the original spatial resolution predicted by the calibrated random forest model.

11. The wide-area multi-level soil property fine mapping method based on spatial downscaling according to claim 1, characterized in that: The specific process of using the spline interpolation method that is most suitable for regularly spaced data to perform spatial interpolation on the soil residuals of the original resolution in step 2) is as follows: 2.1) Determine the interpolation area Calculate the grid size of the high-resolution image, i.e. the number of rows and columns, according to the target resolution, and ensure that the coordinate system and spatial range of the input image are consistent during the interpolation process; 2.2) Apply spline interpolation Extract the pixel values ​​of the original image and their corresponding spatial coordinates and use them as the input points of the spline interpolation. Use the spline interpolation method to generate the interpolation function. Select the two-dimensional spline interpolation to adapt to the two-dimensional spatial characteristics of the raster image. The two-dimensional spline interpolation is a piecewise quadratic polynomial interpolation method, which is used to construct a series of quadratic polynomials between given data points to ensure smoothness in each interval. The goal of the two-dimensional spline interpolation is to find a set of quadratic polynomials so that the function value is continuous at each data point and the derivatives at these points are also continuous. When there are n+1 data points, the quadratic spline interpolation function S(x) is in the following form: For each interval, x i , x i+1 -, define a quadratic polynomial S i (x): S i (x)=a i (x-x i ) 2 +b i (x-x i )+c i (7); In the above formula (7), a i , b i 、c i is the interpolation point condition S i (x i )=y i , Continuity condition S i (x i+1 )=S i+1 (x i+1 ), derivative continuity condition S′ i (x i+1 ) = S′ i+1 (x i+1 ) and boundary conditions S″0(x0)=0 and S″ n-1 (x n )=0 determined coefficient; 2.3) Generate high-resolution images The interpolation results are converted into raster images and then back into raster image format. It is necessary to ensure that the new image maintains the same coordinate system and spatial reference as the target resolution image to obtain the soil property residual map of the target resolution.

Images