Soil environment quality assessment method based on big data

Abstract

The invention relates to the technical field of environment monitoring and big data analysis, and discloses a soil environment quality evaluation method based on big data. The method comprises the following steps: fusing measured data of soil heavy metals with multi-source heterogeneous data such as remote sensing, terrain, land utilization, hydrogeology and the like; a space-attribute joint constraint depth generation network is constructed, an encoder of the space-attribute joint constraint depth generation network extracts multi-scale features, an attribute physical constraint module is embedded into a convection-dispersion equation to model a heavy metal migration mechanism, and a decoder outputs a high-fidelity concentration field; and dividing four pollution risk grades according to national standards, and generating an evaluation diagram. Through data driving and physical mechanism collaborative modeling, the pollution boundary depicting precision is improved, the prediction error under the low sampling density is reduced compared with a traditional method, the recall rate of a high-risk area is increased, and reliable support is provided for precise soil treatment.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of environmental monitoring and big data analysis, and particularly relates to a soil environmental quality evaluation method based on big data. BACKGROUND

[0002] With the continuous improvement of ecological environmental protection requirements, soil environmental quality evaluation has become a key link for pollution prevention and control and sustainable use of land. Among numerous pollutants, heavy metals have become the focus of soil environmental risk control due to their high toxicity, non-degradability and bioaccumulation.

[0003] The current mainstream evaluation method usually relies on field sampling combined with spatial interpolation technology to deduce the regional pollution distribution pattern through limited monitoring point data. Such methods perform well in data-intensive areas, but in wide-scale or complex terrain conditions, due to the constraints of sampling cost and manpower and material resources, only sparse and non-uniformly distributed sample points can be obtained, resulting in uncertainty in the evaluation results.

[0004] Spatial interpolation models based on geographic information systems (GIS) are widely used in spatial prediction of soil heavy metal pollution, such as inverse distance weighting (IDW), Kriging and its variants. These methods assume that the pollutants are continuous and isotropic in space, and their interpolation accuracy is highly dependent on the sampling point density and spatial representativeness. However, in actual soil environment, heavy metal migration and enrichment are affected by multiple factors, including parent material type, land use mode, hydrological conditions and human activities, making the pollution distribution show strong local heterogeneity and boundary mutation characteristics.

[0005] The existing technology has difficulty in effectively capturing such non-stationary spatial variation rules in the case of sparse sampling, especially at the junction of administrative divisions, the transition zone between farmland and industrial areas, or the junction of different types of soil units, often overestimating the pollution level in low-risk areas or underestimating the pollution level in high-risk hotspots. This not only misleads environmental management decisions, but also may cause resource mismatch or risk omission. Although some studies attempt to introduce auxiliary variables (such as remote sensing images, digital elevation models) to improve interpolation accuracy, they lack systematic integration of soil background properties, especially the topological structure and spatial semantic relationship contained in the soil type layer, and fail to fully utilize the constraints and guidance of existing geographic information under sparse data conditions. Therefore, there is an urgent need for a soil environmental quality evaluation method that can deeply integrate soil type topological relationships and adaptively optimize interpolation weights to break through the bottleneck of sparse sampling and achieve high-precision and high-robustness pollution risk spatial characterization. SUMMARY

[0006] This invention provides a big data-based method for assessing soil environmental quality, aiming to solve the technical problem in soil heavy metal pollution assessment where sparse sampling points lead to errors in traditional interpolation methods at pollution boundary regions, resulting in overestimation of low-risk areas or underestimation of high-risk areas. This invention achieves high-fidelity reconstruction of the soil heavy metal concentration field and accurate characterization of risk boundaries by fusing multi-source heterogeneous environmental data, constructing a spatial-attribute joint constraint model, and introducing a deep generative network guided by physical mechanisms.

[0007] According to one aspect of the present invention, a method for assessing soil environmental quality based on big data is provided, comprising: acquiring measured sampling data of soil heavy metals, remote sensing image data, topographic data, land use type data and hydrogeological data within a target area; and performing outlier removal and unit unification processing on the measured sampling data of soil heavy metals to form a standardized sampling point set. Vegetation cover index, surface temperature, and soil moisture information are extracted from the remote sensing image data; slope, aspect, curvature, and runoff accumulation are calculated based on the topographic data; and a land use coding matrix is ​​generated based on the land use type data. Based on the aforementioned hydrogeological data, a groundwater flow direction and permeability coefficient distribution map is constructed. The standardized sampling point set, vegetation cover index, surface temperature, soil moisture, slope, aspect, curvature, runoff accumulation, land use coding matrix, groundwater flow direction and permeability coefficient distribution map are used as input features and input into a pre-trained spatial-attribute joint constraint deep generation network. The spatial-attribute joint constraint deep generation network extracts multi-scale spatial context features through the encoder and combines the attribute physical constraint module to perform mechanism modeling of the heavy metal migration and diffusion process, generating potential feature representations that satisfy the laws of mass conservation and diffusion dynamics. Based on the latent feature representation, the decoder outputs the predicted heavy metal concentration for each spatial cell within the target region; Based on the predicted heavy metal concentrations, pollution risk levels are classified according to national soil environmental quality standards, and a soil environmental quality assessment map is generated.

[0008] Furthermore, the measured soil heavy metal sampling data includes the geographical coordinates of the sampling point, the sampling depth, and the concentration values ​​of five heavy metal elements: cadmium, lead, arsenic, mercury, and chromium. The outlier removal method is three times the standard deviation, and the concentration value sequence of each heavy metal element is performed independently. The unit standardization process converts all concentration values ​​into milligrams per kilogram.

[0009] Furthermore, the remote sensing image data originates from the Gaofen series satellites or the ZY-3 satellite, with a spatial resolution of 2-30 meters; the vegetation cover index is calculated using the normalized difference vegetation index algorithm, with the following formula: ; The surface temperature was obtained by inversion using a single-window algorithm; The soil moisture content was derived using a microwave remote sensing backscattering coefficient and dielectric constant model.

[0010] Furthermore, the topographic data is derived from a digital elevation model, with a spatial resolution of no less than 5 meters; the slope, aspect, and curvature are calculated using the first and second partial derivatives of the digital elevation model; the cumulative flow is calculated using the D8 flow direction algorithm, representing the number of upstream pixels flowing into the pixel per unit area.

[0011] Furthermore, the land use type data comes from the national land use status remote sensing monitoring database, and the classification system includes six categories: cultivated land, forest land, grassland, water area, construction land and unused land; the land use coding matrix is ​​an integer matrix aligned with the spatial grid of the target area, and each pixel is assigned a unique identifier code for its respective land use type.

[0012] Furthermore, the hydrogeological data includes aquifer thickness, porosity, permeability coefficient, and groundwater depth; the groundwater flow direction is determined by Darcy's law combined with regional water level contour maps; and the permeability coefficient distribution map is generated by spatially displaying the measured values ​​of borehole permeability coefficient using the Kriging interpolation method.

[0013] Furthermore, the spatial-attribute joint constraint deep generation network consists of three parts: an encoder, an attribute physical constraint module, and a decoder. The encoder adopts a U-shaped convolutional neural network structure, which includes four downsampling stages and four upsampling stages. Each downsampling stage includes two convolutional layers with a kernel size of 3×3 and a stride of 1, and a max pooling layer with a stride of 2. Each upsampling stage includes a bilinear interpolation upsampling layer and two convolutional layers with a kernel size of 3×3 and a stride of 1. The skip connection concatenates the downsampled feature map and the upsampled feature map at the corresponding scale along the channel dimension; the attribute physical constraint module is embedded after the bottleneck layer of the encoder, receives the high-dimensional feature map output by the encoder and the hydrogeological and topographic features in the aforementioned multi-source environmental features, and simulates the convection-dispersion process of heavy metals in the soil medium through a differentiable finite difference operator, the governing equation of which is: ;in, This represents the concentration of pollutants in groundwater. For time, The groundwater velocity vector. For concentration gradient, For Hamiltonian operators, The dispersion coefficient is a second-order tensor; The dispersion coefficient tensor is determined by the permeability coefficient, porosity, and flow velocity; the module outputs a corrected feature map that satisfies physical laws; the decoder receives the corrected feature map, upsamples it step by step to restore the original spatial resolution, and finally outputs a single-channel heavy metal concentration prediction map.

[0014] Furthermore, the training process of the spatial-attribute joint constraint deep generative network includes: constructing a training sample set, where each sample contains a multi-source environmental feature tensor and the corresponding measured concentration map of heavy metals; The measured heavy metal concentration map is generated by nearest neighbor interpolation of a standardized set of sampling points, retaining only the valid values ​​of the sampling point locations and marking the other locations as missing. The loss function consists of three terms: the first term is the mean square error between the predicted and measured values ​​at the sampling point; the second term is the total variation regularization term of the predicted concentration field, used to suppress unnecessary oscillations; and the third term is the physical residual term, which is the sum of squared residuals calculated on the predicted concentration field by the convection-diffusion equation based on the attribute physical constraint module. The network parameters are optimized by stochastic gradient descent, with an initial learning rate of 0.001, and decay using a cosine annealing strategy.

[0015] Furthermore, the pollution risk levels are divided into four categories: no risk, low risk, medium risk, and high risk. The classification thresholds are based on the risk screening and control values ​​for cadmium, lead, arsenic, mercury, and chromium as specified in the "Soil Environmental Quality Standard for Agricultural Land Soil Pollution Risk Control." For multi-element complex pollution situations, the Nemerow Comprehensive Pollution Index method is used to calculate the comprehensive pollution index, with the formula: ;in The highest single pollution index, It is the arithmetic mean of the squares of each individual pollution index; the comprehensive pollution index value range is mapped to the four risk levels.

[0016] This invention provides a soil environmental quality assessment system based on big data, comprising: a multi-source data acquisition module for acquiring measured sampling data of soil heavy metals, remote sensing image data, topographic data, land use type data, and hydrogeological data within a target area; and a data preprocessing module for removing outliers and standardizing units in the measured sampling data of soil heavy metals to form a standardized sampling point set, and for spatially registering and resampling other multi-source data to give it a unified spatial grid. The feature extraction module is used to extract vegetation cover index, surface temperature, soil moisture, slope, aspect, curvature, runoff accumulation, land use coding matrix, groundwater flow direction, and permeability coefficient distribution map from preprocessed multi-source data; the deep generation and evaluation module is used to input the extracted features into a pre-trained spatial-attribute joint constraint deep generation network to generate predicted heavy metal concentration values ​​for each spatial unit within the target area; the risk level classification module is used to classify pollution risk levels based on the predicted heavy metal concentration values ​​according to national soil environmental quality standards and generate a soil environmental quality assessment map.

[0017] The multi-source data acquisition module connects to the National Earth System Science Data Center, the National Satellite Remote Sensing Application Center of the Ministry of Natural Resources, and the Soil Environment Monitoring Database of the Ministry of Ecology and Environment through a Geographic Information System interface; the data preprocessing module uses bilinear interpolation to resample remote sensing images and topographic data, and the nearest neighbor method to resample land use type data, ensuring that the spatial resolution of all data layers is uniformly 10 meters; the feature extraction module has a built-in standardized algorithm library that supports batch automated calculation; the depth generation and evaluation module is deployed on a distributed computing cluster and supports parallel inference; the risk level classification module has a built-in national standard threshold library that supports dynamic updates.

[0018] Compared with the prior art, the beneficial effects of the present invention are as follows: 1. This invention abandons the traditional evaluation paradigm that relies solely on sparse sampling points for pure mathematical interpolation. Instead, it constructs an evaluation model driven by both data and physical mechanisms by deeply integrating multi-dimensional environmental covariates such as remote sensing, topography, land use, and hydrogeology.

[0019] 2. The spatial-attribute joint constraint deep generative network not only learns the spatial autocorrelation of heavy metal concentration, but also forces the model output to conform to the physical laws of heavy metal migration and diffusion in soil by embedding differentiable convection-diffusion equations. This effectively suppresses oversmoothing or oscillation in sparsely sampled areas, especially in pollution boundary areas.

[0020] 3. This method improves the identification accuracy of high-risk pollution core areas and avoids misjudging low-risk areas, providing reliable technical support for precise soil pollution control and risk management. Experiments show that, under the condition of sampling density of less than 5 points per 100 square kilometers, the root mean square error of the prediction in the pollution boundary area of ​​this invention is lower than that of the ordinary Kriging interpolation method and lower than that of the inverse distance weighting method, and the recall rate of high-risk areas is improved. Attached Figure Description

[0021] Figure 1 This is a schematic diagram of the overall technical solution architecture of the soil environmental quality assessment method based on big data proposed in this invention; Figure 2 This is a schematic diagram of the core principle framework of the spatial-attribute joint constraint deep generation network in this invention; Figure 3 This is a logical flowchart of the multi-source environmental data preprocessing and feature construction in this invention; Figure 4 This is a schematic diagram of the mechanism modeling framework for embedding the attribute physical constraint module into the convection-dispersion equation in this invention; Figure 5 This is a logical flowchart of the heavy metal concentration field reconstruction and pollution risk level classification in this invention. Figure 6 This is a schematic diagram of the multi-level interaction relationship and data flow between the various functional modules of the terminal system and the external data source in this invention. Detailed Implementation

[0022] Please refer to Figures 1 to 6 This invention provides a soil environmental quality assessment method based on big data. Its core is to solve the problem of prediction bias in the pollution boundary area caused by the sparse sampling points of soil heavy metals in traditional spatial interpolation methods.

[0023] Traditional methods such as Kriging interpolation or inverse distance weighting rely solely on the spatial distribution of a limited number of sampling points for mathematical fitting. They neglect the physical mechanisms of migration and diffusion of heavy metals in soil media, which are influenced by multiple environmental factors such as topography, hydrology, land use, and vegetation cover. As a result, under low sampling density conditions, it is easy to overestimate the pollution level in low-risk areas or underestimate the actual degree of harm in high-risk areas.

[0024] To overcome the above-mentioned shortcomings, this invention constructs a deep generative network architecture that integrates multi-source heterogeneous environmental data and physical mechanism guidance. Through spatial-attribute joint constraint modeling, it achieves high-fidelity reconstruction of soil heavy metal concentration field and accurate characterization of pollution risk boundary.

[0025] The method first acquires soil heavy metal sampling data, remote sensing image data, topographic data, land use type data, and hydrogeological data within the target area. The soil heavy metal sampling data includes the geographic coordinates and sampling depth of each sampling point, as well as the concentration values ​​of five heavy metal elements: cadmium, lead, arsenic, mercury, and chromium.

[0026] All concentration values ​​are recorded in milligrams per kilogram. For each heavy metal element concentration sequence, outlier removal is performed independently using the three-standard-deviation method: the mean and standard deviation of the element are calculated across all sampling points, and values ​​exceeding the mean plus or minus three standard deviations are identified as outliers and removed.

[0027] The remaining valid data were then standardized in terms of units to ensure that all concentration values ​​were in the same dimension system, thus forming a standardized set of sampling points.

[0028] The remote sensing image data comes from the Gaofen series satellites or the ZYZ-3 satellite, with a spatial resolution between 2 and 30 meters. Based on this image data, three key surface parameters are extracted: vegetation cover index, surface temperature, and soil moisture.

[0029] The vegetation cover index is calculated using the normalized difference vegetation index algorithm, and its mathematical expression is: ; Surface temperature was obtained through a single-window algorithm, which estimates the true surface temperature using thermal infrared radiance, atmospheric transmittance, and surface emissivity. Soil moisture was derived using a microwave remote sensing backscattering coefficient and dielectric constant model. This model establishes a nonlinear relationship between soil dielectric constant and volumetric water content, and combines radar backscattering intensity to infer the surface soil water content.

[0030] The topographic data is derived from a digital elevation model with a spatial resolution of at least 5 meters. Based on this model, four topographic features are calculated: slope, aspect, curvature, and runoff accumulation.

[0031] Slope represents the degree of inclination of a surface unit, synthesized from the first-order partial derivatives of the digital elevation model in the east-west and north-south directions; aspect represents the direction of maximum descent, obtained by calculating the ratio of the first-order partial derivatives using the arctangent function; curvature reflects the concavity and convexity of the surface morphology, composed of the combination of second-order partial derivatives. The cumulative runoff volume is calculated using the D8 flow direction algorithm, which assumes that the water flows along the steepest slope among the eight neighborhood directions and accumulates the number of pixels that flow into that pixel from upstream, in order to characterize the surface runoff accumulation capacity.

[0032] Land use type data is sourced from the National Land Use Status Remote Sensing Monitoring Database, whose classification system includes six categories: cultivated land, forest land, grassland, water area, construction land, and unused land. To facilitate subsequent model input, this classification data is converted into an integer matrix aligned with the spatial grid of the target area, i.e., a land use coding matrix. Each spatial cell is assigned a unique identifier corresponding to its land use type; for example, cultivated land is coded as 1, forest land as 2, and so on.

[0033] Hydrogeological data includes aquifer thickness, porosity, permeability coefficient, and groundwater depth. Groundwater flow direction is determined using Darcy's law combined with regional water level contour maps: first, water level contour maps are drawn based on borehole water level data, and the hydraulic gradient is calculated; then, the actual groundwater flow direction vector is derived using Darcy's law. The permeability coefficient distribution map is generated spatially from borehole measured permeability coefficient values ​​using Kriging interpolation. This method considers spatial autocorrelation and generates a continuous permeability coefficient field.

[0034] After acquiring the multi-source data described above, data preprocessing steps were performed to ensure that all data layers had a unified spatial grid structure. Remote sensing imagery and topographic data were resampled to 10-meter resolution using bilinear interpolation; land use type data, being categorical variables, were resampled using the nearest neighbor method to avoid category confusion; hydrogeological data were also resampled to the same grid. All data were spatially registered to form aligned multi-channel feature tensors.

[0035] Subsequently, specific feature variables are extracted from the preprocessed multi-source data: vegetation cover index, surface temperature, soil moisture, slope, aspect, curvature, runoff accumulation, land use coding matrix, groundwater flow direction vector field, and permeability coefficient distribution map. These features together constitute the input of the deep generative network.

[0036] The deep generative network is a pre-trained spatial-attribute joint constraint deep generative network, which consists of three parts: an encoder, an attribute physical constraint module, and a decoder. The encoder adopts a U-shaped convolutional neural network structure, including four downsampling stages and four upsampling stages.

[0037] Each downsampling stage consists of two convolutional layers with a kernel size of 3×3 and a stride of 1, followed by a max-pooling layer with a stride of 2, used to progressively extract multi-scale spatial context features. Each upsampling stage consists of a bilinear interpolation upsampling layer and two convolutional layers with a kernel size of 3×3 and a stride of 1. Skip connections concatenate the downsampling feature maps and upsampling feature maps of the corresponding scale along the channel dimension, preserving fine-grained spatial information.

[0038] The attribute physical constraint module is embedded after the bottleneck layer of the encoder. This module receives the high-dimensional feature map output by the encoder, and simultaneously incorporates groundwater flow direction, permeability coefficient, porosity, slope, and runoff accumulation from hydrogeological and topographic features. Internally, the module constructs a differentiable finite difference operator to simulate the convection-dispersion process of heavy metals in the soil medium. Its governing equations are: ; in, This represents the concentration of pollutants in groundwater. For time, The groundwater velocity vector. For concentration gradient, For Hamiltonian operators, Let be the second-order tensor of the dispersion coefficient; according to Darcy's law Calculations show that Permeability coefficient, For the water head, Porosity; The dispersion coefficient tensor has its principal axis aligned with the flow velocity direction, and its magnitude is determined by the permeability coefficient, porosity, and flow velocity. This equation is implemented on a discrete grid using a central difference scheme to calculate whether the predicted concentration field satisfies mass conservation and diffusion kinetics. If not, a physical residual signal is generated to correct the feature representation. The module outputs a corrected feature map that satisfies the physical laws.

[0039] The decoder receives the corrected feature map, upsamples it step by step to restore it to the original 10-meter spatial resolution, and finally outputs a single-channel heavy metal concentration prediction map, covering every spatial unit in the target area.

[0040] The training process of the spatial-attribute joint constraint deep generative network is as follows: A training sample set is constructed, where each sample contains a multi-source environmental feature tensor and a corresponding measured heavy metal concentration map. The measured concentration map is generated by nearest-neighbor interpolation of a standardized set of sampling points, retaining only valid concentration values ​​at the sampling point locations, and marking the remaining locations as missing values.

[0041] The loss function consists of three terms: The first term is the mean square error between the predicted and measured values ​​at the sampling point location. , denoted as: ; in This represents the number of valid sampling points. For the first Predicted concentration values ​​at each sampling point; For the first The measured concentration values ​​at each sampling point.

[0042] The second term is the regularization term for the total variation of the predicted concentration field. : ; Used to suppress unnecessary oscillations and enhance spatial smoothness; To predict the concentration field in The absolute value of the gradient in the direction; To predict the concentration field in The absolute value of the gradient in the direction.

[0043] The third term is the physical residual. : ; To predict concentration; to measure the deviation of the predicted field from the convection-diffusion equation. The total loss function is: ,in and The weights are set to 0.01 and 0.1 respectively. The network parameters are optimized using stochastic gradient descent with an initial learning rate of 0.001, which is then decayed to 0.0001 using cosine annealing.

[0044] After completing the concentration field reconstruction, the risk level classification stage begins. Based on the risk screening and control values ​​for cadmium, lead, arsenic, mercury, and chromium as specified in the "Soil Environmental Quality Standard for Agricultural Land Soil Pollution Risk Control," a separate pollution index is calculated for each heavy metal. ,in To predict concentration, This is the risk screening value for the corresponding element.

[0045] For multi-element complex pollution situations, the Nemerow comprehensive pollution index method is used to calculate the comprehensive pollution index. Its formula is: ; in The highest single pollution index, This is the arithmetic mean of the squares of each individual pollution index. According to... The numerical range is used to classify pollution risks into four categories: when At that time, it was determined to be risk-free; when At that time, it was determined to be low risk; when At that time, it was determined to be a medium-risk area; when At that time, it was determined to be high risk.

[0046] Based on this, a soil environmental quality assessment map is generated, with each spatial unit in the map labeled with its corresponding risk level.

[0047] Throughout the methodology, data flow strictly adheres to the principle of spatiotemporal consistency. All multi-source data undergoes spatial registration and temporal alignment before being input into the network, ensuring that feature variables are physically comparable. During network inference, the attribute physics constraint module verifies the physical rationality of the concentration field in real time, preventing the model from generating non-physical interpretations due to data sparsity. The risk classification stage employs a national standard threshold library, supporting dynamic updates to adapt to the latest environmental regulations.

[0048] The proposed method maintains high accuracy even under extremely sparse sampling conditions, with a sampling density of less than 5 points per 100 square kilometers. Experimental verification shows that, compared to the ordinary Kriging interpolation method, the proposed method reduces the root mean square error in the prediction of polluted boundary regions; compared to the inverse distance weighting method, the error is also reduced. Simultaneously, the recall rate in high-risk areas is improved, reducing false negatives and false negatives, and providing a reliable decision-making basis for precise soil pollution control.

[0049] As a system implementation of the present invention, the big data-based soil environmental quality assessment system includes a multi-source data acquisition module, a data preprocessing module, a feature extraction module, a deep generation assessment module, and a risk level classification module.

[0050] The multi-source data acquisition module connects to the National Earth System Science Data Center, the National Satellite Remote Sensing Application Center of the Ministry of Natural Resources, and the Soil Environmental Monitoring Database of the Ministry of Ecology and Environment via a Geographic Information System interface, automatically downloading the required data for the target area. The data preprocessing module performs spatial registration, resampling, and format conversion, unifying all data layers to a 10-meter grid.

[0051] The feature extraction module incorporates a standardized algorithm library to batch calculate derived variables such as vegetation cover index, slope, and runoff accumulation. The deep generation and assessment module is deployed on a distributed computing cluster, supporting parallel inference across multiple regions, with a single inference cycle taking less than 15 minutes. The risk level classification module integrates a national standard threshold library, automatically generating a four-color risk assessment map and supporting export to a GIS-compatible format.

[0052] During system operation, each module executes sequentially: the multi-source data acquisition module starts and simultaneously retrieves five types of data; the data preprocessing module performs geometric correction and resolution unification on the raw data; the feature extraction module generates an eleven-dimensional feature tensor; the deep learning and evaluation module loads a pre-trained model, inputs the feature tensor, and outputs a heavy metal concentration prediction map; the risk level classification module reads the prediction map and threshold library, completes risk classification, and generates the final evaluation product. The entire process requires no manual intervention and has end-to-end automation capabilities.

[0053] In summary, this invention, through a deep learning architecture that deeply integrates multi-source environmental covariates and physical mechanisms, breaks through the performance bottleneck of traditional interpolation methods under sparse sampling conditions, and achieves high-precision reconstruction of soil heavy metal pollution fields and accurate identification of risk boundaries, providing a scientific and reliable technical tool for ecological and environmental management.

Claims

1. A soil environmental quality assessment method based on big data, characterized in that, include: Acquire measured sampling data of heavy metals in the soil, remote sensing image data, topographic data, land use type data, and hydrogeological data within the target area; The measured soil heavy metal sampling data were processed to remove outliers and standardize units to form a standardized sampling point set; Based on the remote sensing image data, vegetation cover index, surface temperature and soil moisture information are extracted; Calculate the slope, aspect, curvature, and cumulative runoff based on the topographic data; A land use coding matrix is ​​generated based on the land use type data; Based on the aforementioned hydrogeological data, a map showing the direction of groundwater flow and the distribution of permeability coefficients was constructed. The standardized sampling point set, vegetation cover index, surface temperature, soil moisture, slope, aspect, curvature, runoff accumulation, land use coding matrix, groundwater flow direction and permeability coefficient distribution map are used as input features and input into a pre-trained spatial-attribute joint constraint deep generation network. The spatial-attribute joint constraint deep generation network extracts multi-scale spatial context features through the encoder and combines the attribute physical constraint module to perform mechanism modeling of the heavy metal migration and diffusion process, generating potential feature representations that satisfy the laws of mass conservation and diffusion dynamics. Based on the latent feature representation, the decoder outputs the predicted heavy metal concentration for each spatial cell within the target region; Based on the predicted heavy metal concentrations, pollution risk levels are classified according to national soil environmental quality standards, and a soil environmental quality assessment map is generated.

2. The soil environmental quality assessment method based on big data according to claim 1, characterized in that, The measured soil heavy metal sampling data were processed to remove outliers and standardize units, forming a standardized sampling point set, including: The measured soil heavy metal sampling data includes the geographical coordinates of the sampling point, the sampling depth, and the concentration values ​​of five heavy metal elements: cadmium, lead, arsenic, mercury, and chromium. Outlier removal was performed independently on the concentration value sequence of each heavy metal element using the three-standard-deviation method. All concentration values ​​were uniformly converted to milligrams per kilogram to form a standardized sampling point set.

3. The soil environmental quality assessment method based on big data according to claim 2, characterized in that, Based on the remote sensing image data, vegetation cover index, surface temperature, and soil moisture information are extracted, including: The remote sensing image data comes from the Gaofen series satellites or the Ziyuan-3 satellite, with a spatial resolution of 2-30 meters. The vegetation cover index is calculated using the normalized difference vegetation index algorithm, and the formula is as follows: ； The surface temperature was obtained by inversion using a single-window algorithm; The soil moisture content was derived using a microwave remote sensing backscattering coefficient and dielectric constant model.

4. The soil environmental quality assessment method based on big data according to claim 3, characterized in that, Based on the aforementioned topographic data, the slope, aspect, curvature, and cumulative runoff volume are calculated, including: The terrain and geomorphology data are derived from a digital elevation model, and its spatial resolution is no less than five meters. The slope, aspect, and curvature are calculated using the first and second partial derivatives of the digital elevation model. The cumulative flow is calculated using the D8 flow direction algorithm and represents the number of upstream pixels flowing into the pixel per unit area.

5. The soil environmental quality assessment method based on big data according to claim 4, characterized in that, Based on the land use type data, a land use coding matrix is ​​generated, including: The land use type data comes from the National Land Use Status Remote Sensing Monitoring Database. The classification system includes six categories: cultivated land, forest land, grassland, water area, construction land, and unused land. The land use coding matrix is ​​an integer matrix aligned with the spatial grid of the target area, and each cell is assigned a unique identifier for its respective land use type.

6. The soil environmental quality assessment method based on big data according to claim 5, characterized in that, Based on the aforementioned hydrogeological data, a map of groundwater flow direction and permeability coefficient distribution is constructed, including: The hydrogeological data include aquifer thickness, porosity, permeability coefficient, and groundwater depth; The direction of groundwater flow was determined by Darcy's law in conjunction with a regional water level contour map. The permeability coefficient distribution map is generated by spatially displaying the measured values ​​of borehole permeability coefficient using the Kriging interpolation method.

7. The soil environmental quality assessment method based on big data according to claim 6, characterized in that, The spatial-attribute joint constraint deep generation network extracts multi-scale spatial context features through an encoder and combines them with an attribute physical constraint module to perform mechanistic modeling of the heavy metal migration and diffusion process, generating potential feature representations that satisfy mass conservation and diffusion dynamics, including: The encoder adopts a U-shaped convolutional neural network structure, which includes four downsampling stages and four upsampling stages. Each downsampling stage includes two convolutional layers with a kernel size of 3×3 and a stride of 1, and a max pooling layer with a stride of 2. Each upsampling stage includes a bilinear interpolation upsampling layer and two convolutional layers with a kernel size of 3×3 and a stride of 1. Skip connections concatenate the downsampling feature maps and upsampling feature maps of the corresponding scale along the channel dimension. The attribute physical constraint module is embedded after the bottleneck layer of the encoder. It receives the high-dimensional feature map and hydrogeological and topographic features output by the encoder, and simulates the convection-dispersion process of heavy metals in the soil medium through a differentiable finite difference operator. The governing equation is: ;in, This represents the concentration of pollutants in groundwater. For time, The groundwater velocity vector. For concentration gradient, For Hamiltonian operators, The dispersion coefficient is a second-order tensor; the dispersion coefficient tensor is determined by the permeability coefficient, porosity, and flow velocity; the module outputs a corrected feature map that satisfies physical laws.

8. The soil environmental quality assessment method based on big data according to claim 7, characterized in that, The training process of the spatial-attribute joint constraint deep generative network includes: A training sample set is constructed, and each sample contains a multi-source environmental feature tensor and a corresponding measured heavy metal concentration map. The measured heavy metal concentration map is generated by nearest neighbor interpolation of a standardized sampling point set, retaining only the valid values ​​of the sampling point positions and marking the other positions as missing. The loss function consists of three terms: the first term is the mean square error between the predicted and measured values ​​at the sampling point location; the second term is the total variation regularization term of the predicted concentration field; and the third term is the physical residual term, which is the sum of squared residuals calculated on the predicted concentration field by the convection-diffusion equation based on the attribute physical constraint module. The network parameters were optimized using stochastic gradient descent with an initial learning rate of 0.001, which was decayed using a cosine annealing strategy.

9. The soil environmental quality assessment method based on big data according to claim 8, characterized in that, Based on the predicted heavy metal concentrations, pollution risk levels are classified according to national soil environmental quality standards, and a soil environmental quality assessment map is generated, including: The pollution risk levels are divided into four categories: no risk, low risk, medium risk, and high risk. The threshold values ​​are based on the risk screening and control values ​​for cadmium, lead, arsenic, mercury, and chromium specified in the "Soil Environmental Quality Agricultural Land Soil Pollution Risk Control Standard". For multi-element complex pollution situations, the Nemerow Comprehensive Pollution Index method is used to calculate the comprehensive pollution index, and the formula is as follows: ;in The highest single pollution index, It is the arithmetic mean of the squares of each individual pollution index; the comprehensive pollution index value range is mapped to the four risk levels.

10. The soil environmental quality assessment method based on big data according to claim 9, characterized in that, All multi-source data undergo spatial registration and resampling before being input into the spatial-attribute joint constraint deep generation network, so that it has a uniform spatial grid. The remote sensing images and terrain data were resampled to 10-meter resolution using bilinear interpolation. The land use type data were resampled using the nearest neighbor method. The hydrogeological data are resampled to the same grid to form aligned multi-channel feature tensors.

Images