Non-point source phosphorus pollution prediction method based on small watershed in mountainous area

Abstract

The invention relates to the field of non-point source phosphorus pollution prediction, in particular to a non-point source phosphorus pollution prediction method based on a small watershed in a mountainous area. The scheme comprises the following steps: collecting basic geographic data, meteorological data and hydrology and water quality data of a small watershed of a target mountainous area, and performing watershed division and water system generation by using a digital elevation model; based on unmanned aerial vehicle multispectral soil key parameter inversion, a rainfall erosion force factor, a soil erodibility factor, a slope length factor, a vegetation coverage and management factor and a water and soil conservation measure factor are calculated by utilizing unmanned aerial vehicle remote sensing data, and the factors are multiplied to generate a watershed high-resolution soil erosion intensity distribution diagram. Meanwhile, a soil total phosphorus distribution diagram is generated, and an improved SWAT model is constructed; carrying out pollution load spatialization of Kriging interpolation on the basis of an improved SWAT model of completion calibration and verification; and finally, predicting the non-point source phosphorus pollution load in a future set scene. The method is suitable for non-point source phosphorus pollution prediction of the small watershed in the mountainous area.

Description

Technical Field

[0001] This invention relates to the field of non-point source phosphorus pollution prediction, and specifically to a method for predicting non-point source phosphorus pollution based on small watersheds in mountainous areas. Background Technology

[0002] Non-point source pollution is one of the most important issues in current water environment protection, especially phosphorus pollution, which is considered a major driver of eutrophication in lakes and rivers. Unlike point source pollution, non-point source pollution is characterized by its wide distribution, complex sources, and large spatiotemporal variability, and is particularly prominent in small watersheds in mountainous areas. The natural conditions of small watersheds in mountainous areas (such as complex topography, high proportion of sloping farmland, and diverse land use types) lead to extremely complex dynamic responses in the rainfall-runoff-erosion-migration process, making it difficult to predict the migration of pollutants, especially phosphorus.

[0003] In traditional research, non-point source pollution prediction mainly relies on the following types of methods:

[0004] Mechanistic modeling: The representative model is the SWAT (Soil and Water Assessment Tool) model. The SWAT model is a comprehensive simulation model based on watershed physical processes, capable of quantitatively simulating the entire process from rainfall, runoff generation, erosion, sediment transport to pollutant migration. This model typically includes data preparation, sub-watershed delineation, hydrological response unit generation, hydrological cycle and runoff generation calculation, sediment yield calculation (based on the MUSLE equation), and a pollutant migration module. However, the SWAT model is extremely sensitive to the resolution and accuracy of the input data. In small mountainous watersheds, due to steep terrain, scarce data, and strong spatial heterogeneity, model parameter calibration is difficult, and simulation accuracy is often insufficient.

[0005] Statistical empirical modeling: This method establishes the relationship between pollutant output and environmental factors based on measured data through algorithms such as regression analysis, neural networks, or support vector machines. While computationally simple, this method lacks a physical basis, cannot explain pollution migration mechanisms, and has poor predictive ability in areas with scarce data.

[0006] On-site monitoring method: High-precision information is obtained by setting up monitoring points and collecting high-frequency data. However, due to inconvenient transportation and strong spatial heterogeneity in mountainous watersheds, long-term monitoring is costly, has a long cycle, and limited coverage.

[0007] Remote sensing estimation: Satellite remote sensing (such as Landsat and Sentinel-2) can be used for land use and vegetation estimation, but it is limited by resolution and cloud cover, making it difficult to support high-precision modeling of small-scale watersheds. Summary of the Invention

[0008] The purpose of this invention is to overcome the shortcomings of the prior art and provide a non-point source phosphorus pollution prediction method based on small watersheds in mountainous areas, thereby improving the prediction accuracy.

[0009] The present invention achieves the above objectives by adopting the following technical solution: The present invention provides a method for predicting non-point source phosphorus pollution based on small watersheds in mountainous areas, comprising:

[0010] S1. Acquire watershed data and perform preprocessing;

[0011] Collect basic geographic data, meteorological data, and hydrological and water quality data for the target mountainous small watershed, and use digital elevation models to divide the watershed and generate the river system;

[0012] S2. Inversion of key soil parameters based on UAV multispectral data;

[0013] Soil samples were collected in a gridded manner in designated plots within the watershed, and the total phosphorus content and soil organic matter content were measured simultaneously.

[0014] Simultaneously with soil sampling, drones equipped with multispectral and lidar sensors were used to fly over the entire watershed, acquiring high-resolution multispectral images and high-precision digital surface models.

[0015] Radiometric calibration, atmospheric correction, and geometric fine correction are performed on the acquired multispectral images, and vegetation indices are calculated. The lidar data is processed to generate a high-precision digital elevation model.

[0016] Correlation analysis was performed on the measured data of soil sampling points with the corresponding multispectral band values ​​and calculated vegetation indices. Sensitive bands and indices were selected as independent variables. At the same time, the rainfall erosivity factor, soil erodibility factor, slope and slope length factor, vegetation cover and management factor, and soil and water conservation measure factor were calculated using UAV remote sensing data through the general soil loss equation.

[0017] By multiplying the various factors, a high-resolution distribution map of soil erosion intensity in the watershed is generated.

[0018] Correlation analysis was performed on the measured data of soil sampling points with the corresponding multispectral band values ​​and calculated indices to establish an inversion model for total phosphorus in the soil and generate a distribution map of total phosphorus in the soil.

[0019] S3. Construct an improved SWAT model;

[0020] An improved SWAT model was constructed, which directly embeds slope and slope length factors, vegetation cover and management factors into the corresponding calculation modules of the SWAT model, replacing the original estimation method of the SWAT model based on the average value of sub-basins.

[0021] High-resolution soil erosion intensity distribution maps of the watershed were used as additional validation data for SWAT model calibration, and soil total phosphorus distribution maps were used as input to the SWAT model to train the improved SWAT model.

[0022] After training, the hydrological, sediment, and phosphorus migration modules of the improved SWAT model were comprehensively calibrated and validated using measured flow, sediment concentration, and total phosphorus concentration data at the watershed outlet section.

[0023] S4. Based on the improved SWAT model with completion rate determination and verification, spatialize the pollution load using Kriging interpolation;

[0024] S5. Predict the non-point source phosphorus pollution load under future specified scenarios.

[0025] Furthermore, the calculation of rainfall erosivity factors, soil erodibility factors, slope and slope length factors, vegetation cover and management factors, and soil and water conservation measures factors using UAV remote sensing data specifically includes:

[0026] Rainfall erosivity factor: Collect rain gauge data from the watershed and surrounding rain gauges, and calculate the rainfall erosivity factor using the rainfall data;

[0027] Soil erodibility factors are calculated using soil mechanical composition and organic matter data, combined with the EPIC (Erosion-Productivity Impact Calculator) model or the nomograph method, and spatially extended using soil properties retrieved from multispectral data.

[0028] The slope length factor is derived by using raster operations to accurately calculate the slope length and slope of each pixel based on the generated high-precision digital elevation model.

[0029] Vegetation cover and management factors: Using vegetation cover obtained by inversion from multispectral data, a quantitative relationship is established with vegetation cover and management factors, and a high-resolution vegetation cover and management factor map is generated.

[0030] Soil and water conservation measures factors are identified through manual interpretation of UAV orthophotos or deep learning models to determine soil and water conservation measures for terraced fields and contour farming, and corresponding soil and water conservation measures factor values ​​are assigned.

[0031] Furthermore, step S4 specifically includes:

[0032] The improved SWAT model was run to complete the calibration and validation, and output the non-point source phosphorus pollution load for each hydrological response unit within a specified time period.

[0033] The phosphorus pollution load value of each hydrological response unit is assigned to its centroid point to form a discrete pollution load point dataset.

[0034] Using the Kriging interpolation method, spatial interpolation calculations are performed based on the obtained discrete pollution load point dataset to generate a continuous spatial distribution map of non-point source phosphorus pollution load covering the entire watershed.

[0035] Furthermore, step S5 specifically includes:

[0036] Define future projection scenarios, including climate change scenarios and land use change scenarios;

[0037] Prepare driving data for future scenarios. For climate change scenarios, input future climate prediction data. For land use change scenarios, input future land use maps. Update vegetation cover and management factors and soil and water conservation measures factors in the SWAT model in sync based on future land use types.

[0038] Run the improved SWAT model to simulate the hydrological, sediment, and phosphorus migration processes in the watershed over a specified future time period;

[0039] For the phosphorus load of each hydrological response unit output by the improved SWAT model, repeat step S4 to generate a predicted distribution map and trend map of future non-point source phosphorus pollution load.

[0040] The beneficial effects of this invention are as follows:

[0041] This invention utilizes UAV multispectral and lidar data to retrieve high-resolution parameters such as total phosphorus, organic matter, slope, and slope length in the soil, replacing traditional low-resolution input data and improving model accuracy from the source.

[0042] This invention enables rapid acquisition of comprehensive data using drones, reducing the workload of field sampling and laboratory analysis, and achieving high-efficiency, low-cost data updates.

[0043] This invention spatializes the output of the mechanistic model through Kriging interpolation to generate a phosphorus pollution load distribution map, thereby achieving spatial visualization of the migration process.

[0044] This invention embeds high-resolution RUSLE factors (slope length factor, cover management factor, and protection measure factor) into the SWAT model, improves the sediment and particulate phosphorus migration module, and enhances the reliability of the model in complex mountainous areas.

[0045] This invention simulates phosphorus load change trends under different scenarios by inputting future climate and land use data, providing quantitative decision support for ecological restoration and non-point source pollution control. Attached Figure Description

[0046] Figure 1 This is a flowchart of a non-point source phosphorus pollution prediction method based on a small watershed in a mountainous area, provided by an embodiment of the present invention. Detailed Implementation

[0047] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings.

[0048] This invention provides a method for predicting non-point source phosphorus pollution in small watersheds in mountainous areas, such as... Figure 1 As shown, it specifically includes:

[0049] S1. Acquire watershed data and perform preprocessing;

[0050] We collected basic geographic, meteorological, and hydrological data on small watersheds in the target mountainous area, and used digital elevation models to delineate watersheds and generate river systems, providing a foundation for the construction of SWAT models.

[0051] S2. Inversion of key soil parameters based on UAV multispectral data;

[0052] Soil samples were collected in a gridded manner in designated plots within the watershed, and the total phosphorus content and soil organic matter content were measured simultaneously.

[0053] Simultaneously with soil sampling, drones equipped with multispectral and lidar sensors were used to fly over the entire watershed, acquiring high-resolution multispectral images and high-precision digital surface models.

[0054] Radiometric calibration, atmospheric correction, and geometric fine correction are performed on the acquired multispectral images, and vegetation indices are calculated. The lidar data is processed to generate a high-precision digital elevation model.

[0055] Correlation analysis was performed on the measured data of soil sampling points with the corresponding multispectral band values ​​and calculated vegetation indices. Sensitive bands and indices were selected as independent variables. At the same time, the rainfall erosivity factor, soil erodibility factor, slope and slope length factor, vegetation cover and management factor, and soil and water conservation measure factor were calculated using UAV remote sensing data through the general soil loss equation.

[0056] Rainfall erosivity factor: Calculated by collecting data from rain gauge stations in the watershed and surrounding areas, and using rainfall data.

[0057] Soil erodibility factor: Calculated using measured soil mechanical composition and organic matter data, combined with the EPIC model or nomograph method, and spatially extended using soil properties retrieved from multispectral data.

[0058] Slope length factor: Based on the generation of a high-precision digital elevation model, the slope length and slope of each pixel are accurately calculated using raster operations, thereby obtaining the slope length factor map.

[0059] Vegetation cover and management factors: Using vegetation cover obtained from multispectral data inversion, a quantitative relationship is established with vegetation cover and management factors, and a high-resolution vegetation cover and management factor map is generated.

[0060] Soil and water conservation measures factors: Soil and water conservation measures such as terraced fields and contour farming are identified through manual interpretation of UAV orthophotos or deep learning models, and corresponding soil and water conservation measures factor values ​​are assigned.

[0061] Multiply the above factors to generate a high-resolution soil erosion intensity distribution map of the watershed.

[0062] Correlation analysis was performed on the measured data of soil sampling points with the corresponding multispectral band values ​​and calculated indices to establish an inversion model for total phosphorus in the soil and generate a distribution map of total phosphorus in the soil.

[0063] The improvement in this step is that the pixel is changed to an area-weighted direct assignment of hydrological response unit, the slope and slope length factors are calculated by LiDAR micro-topography, the vegetation cover and management factors are calculated by the spectral coverage model, the soil and water conservation measures factors are identified by soil and water conservation measures imagery, and the soil erodibility is estimated by soil properties / spectral data.

[0064] The direct technical effects include: reducing scale mismatch and surface homogenization bias, while preserving micro-topography and cover heterogeneity.

[0065] Verifiable indicators or scenarios: Improved spatial contrast of erosion hotspots; reduced deviation between hydrological response unit results and sample point or area monitoring.

[0066] S3. Construct an improved SWAT model;

[0067] An improved SWAT model was constructed, which directly embeds slope and slope length factors, vegetation cover and management factors into the corresponding calculation modules of the SWAT model, replacing the original estimation method of the SWAT model based on the average value of sub-basins.

[0068] High-resolution soil erosion intensity distribution maps of the watershed were used as additional validation data for SWAT model calibration, and soil total phosphorus distribution maps were used as input to the SWAT model to train the improved SWAT model.

[0069] After training, the hydrological, sediment, and phosphorus migration modules of the improved SWAT model were comprehensively calibrated and validated using measured flow, sediment concentration, and total phosphorus concentration data at the watershed outlet section, ensuring that the model's simulation of water volume, sediment volume, and phosphorus load all met reliable standards.

[0070] The first mechanism improvement in this step is to directly drive the mechanistic model (instead of a quadratic average) with high-resolution rainfall erosivity factors, soil erosivity factors, slope length factors, vegetation cover and management factors, or soil and water conservation measures factors.

[0071] Direct technical effect: The model's response sensitivity to local governance measures (strips, terraces, changes in coverage) is enhanced.

[0072] Verifiable indicators or scenarios: Before and after the implementation of small-scale soil conservation measures, the simulated elasticity of sediment or total phosphorus is significantly enhanced.

[0073] The second improvement in this step involves tightly coupling the RUSLE (Revised Universal Soil Loss Equation) factor field with the SWAT sediment-particle phosphorus equation. The basic principle of RUSLE is to estimate the rate of soil erosion by considering the main driving factors of soil erosion. Its calculation formula is as follows:

[0074] A = R×K×LS×C×P;

[0075] Where A represents soil erosion per unit area, R represents rainfall erosivity factor, K represents soil erodibility factor, LS represents slope and slope length factor, C represents vegetation cover and management factor, and P represents soil and water conservation measures factor.

[0076] Direct technical effect: It closes the "phosphorus follows mud" chain from a mechanistic perspective, avoiding mismatch between mud and sand by only matching water volume or total phosphorus.

[0077] Verifiable indicators or scenarios: The goodness of fit among flow rate, suspended sediment concentration (or sediment flux), and total phosphorus is more balanced, and the situation is more stable during the flood peak or receding period.

[0078] The third mechanism improvement in this step is to introduce a joint calibration of three objectives—flow rate, suspended sediment concentration (or sediment flux), and total phosphorus—with SUFI-2 / GLUE uncertainty constraints.

[0079] Direct technical effect: Improves parameter identifiability and interval reliability.

[0080] Verifiable metrics or scenarios: convergence in the 95% prediction uncertainty band (95PPU); more robust generalization in extrapolation scenarios.

[0081] S4. Based on the improved SWAT model with completion rate determination and verification, spatialize the pollution load using Kriging interpolation;

[0082] The SWAT model with the running rate set outputs the non-point source phosphorus pollution load for each hydrological response unit within a specified time period.

[0083] The phosphorus pollution load value of each hydrological response unit is assigned to its centroid point to form a discrete pollution load point dataset.

[0084] Using the Kriging interpolation method, spatial interpolation calculations are performed based on the obtained discrete pollution load point dataset to generate a continuous spatial distribution map of non-point source phosphorus pollution load covering the entire watershed.

[0085] The mechanism improvement in this step is to use the centroid of the hydrological response unit as the sampling point to perform spatial statistical representation of (co-)kriging interpolation.

[0086] Direct technical effect: Obtain an unbiased, minimum variance continuous risk surface to support threshold classification and governance prioritization.

[0087] Verifiable metrics or scenarios: Hotspot identification is not sensitive to sampling layout or variation function selection, and the partitioning scheme is reproducible.

[0088] S5. Predict the non-point source phosphorus pollution load under future specified scenarios.

[0089] Define future projection scenarios, including climate change scenarios and land use change scenarios;

[0090] Prepare driving data for future scenarios. For climate change scenarios, input future climate prediction data. For land use change scenarios, input future land use maps. Update vegetation cover and management factors and soil and water conservation measures factors in the SWAT model in sync based on future land use types.

[0091] Run the improved SWAT model to simulate the hydrological, sediment, and phosphorus migration processes in the watershed over a specified future time period;

[0092] For the phosphorus load of each hydrological response unit output by the improved SWAT model, repeat step S4 to generate a predicted distribution map and trend map of future non-point source phosphorus pollution load.

[0093] The mechanism improvement in this step is to simultaneously update vegetation cover and management factors or soil and water conservation measures factors and link them with the mechanism model under climate or land use scenarios.

[0094] Direct technical effects: Enables spatiotemporal evolution prediction of climate dual scenarios, supporting the "assessment before investment" approach.

[0095] Verifiable indicators or scenarios: Erosion or total phosphorus elasticity curves and critical thresholds under different typical concentration pathways or shared socio-economic pathways and governance combinations can be quantitatively compared.

[0096] The above description is merely a preferred embodiment of the present invention. It should be understood that the present invention is not limited to the forms disclosed herein and should not be construed as excluding other embodiments. It can be used in various other combinations, modifications, and environments, and can be altered within the scope of the concept described herein through the above teachings or related technologies or knowledge. Modifications and variations made by those skilled in the art that do not depart from the spirit and scope of the present invention should be within the protection scope of the appended claims.

Claims

1. A method for predicting non-point source phosphorus pollution based on small watersheds in mountainous areas, characterized in that, include: S1. Acquire watershed data and perform preprocessing; Collect basic geographic data, meteorological data, and hydrological and water quality data for the target mountainous small watershed, and use digital elevation models to divide the watershed and generate the river system; S2. Inversion of key soil parameters based on UAV multispectral data; Soil samples were collected in a gridded manner in designated plots within the watershed, and the total phosphorus content and soil organic matter content were measured simultaneously. Simultaneously with soil sampling, drones equipped with multispectral and lidar sensors were used to fly over the entire watershed, acquiring high-resolution multispectral images and high-precision digital surface models. Radiometric calibration, atmospheric correction, and geometric fine correction are performed on the acquired multispectral images, and vegetation indices are calculated. The lidar data is processed to generate a high-precision digital elevation model. Correlation analysis was performed on the measured data of soil sampling points with the corresponding multispectral band values ​​and calculated vegetation indices. Sensitive bands and indices were selected as independent variables. At the same time, the rainfall erosivity factor, soil erodibility factor, slope and slope length factor, vegetation cover and management factor, and soil and water conservation measure factor were calculated using UAV remote sensing data through the general soil loss equation. By multiplying the various factors, a high-resolution distribution map of soil erosion intensity in the watershed is generated. Correlation analysis was performed on the measured data of soil sampling points with the corresponding multispectral band values ​​and calculated indices to establish an inversion model for total phosphorus in the soil and generate a distribution map of total phosphorus in the soil. S3. Construct an improved SWAT model; An improved SWAT model was constructed, which directly embeds slope and slope length factors, vegetation cover and management factors into the corresponding calculation modules of the SWAT model, replacing the original estimation method of the SWAT model based on the average value of sub-basins. High-resolution soil erosion intensity distribution maps of the watershed were used as additional validation data for SWAT model calibration, and soil total phosphorus distribution maps were used as input to the SWAT model to train the improved SWAT model. After training, the hydrological, sediment, and phosphorus migration modules of the improved SWAT model were comprehensively calibrated and validated using measured flow, sediment concentration, and total phosphorus concentration data at the watershed outlet section. S4. Based on the improved SWAT model with completion rate determination and verification, spatialize the pollution load using Kriging interpolation; S5. Predict the non-point source phosphorus pollution load under future specified scenarios.

2. The method for predicting non-point source phosphorus pollution based on small watersheds in mountainous areas according to claim 1, characterized in that, The calculation of rainfall erosivity factors, soil erodibility factors, slope and slope length factors, vegetation cover and management factors, and soil and water conservation measures factors using UAV remote sensing data specifically includes: Rainfall erosivity factor: Collect rain gauge data from the watershed and surrounding rain gauges, and calculate the rainfall erosivity factor using the rainfall data; Soil erodibility factors are calculated using soil mechanical composition and organic matter data, combined with the EPIC model or nomograph method, and spatially extended using soil properties retrieved from multispectral data. The slope length factor is derived by using raster operations to accurately calculate the slope length and slope of each pixel based on the generated high-precision digital elevation model. Vegetation cover and management factors: Using vegetation cover obtained by inversion from multispectral data, a quantitative relationship is established with vegetation cover and management factors, and a high-resolution vegetation cover and management factor map is generated. Soil and water conservation measures factors are identified through manual interpretation of UAV orthophotos or deep learning models to determine soil and water conservation measures for terraced fields and contour farming, and corresponding soil and water conservation measures factor values ​​are assigned.

3. The method for predicting non-point source phosphorus pollution based on small watersheds in mountainous areas according to claim 2, characterized in that, Step S4 specifically includes: The improved SWAT model was run to complete the calibration and validation, and output the non-point source phosphorus pollution load for each hydrological response unit within a specified time period. The phosphorus pollution load value of each hydrological response unit is assigned to its centroid point to form a discrete pollution load point dataset. Using the Kriging interpolation method, spatial interpolation calculations are performed based on the obtained discrete pollution load point dataset to generate a continuous spatial distribution map of non-point source phosphorus pollution load covering the entire watershed.

4. The method for predicting non-point source phosphorus pollution based on small watersheds in mountainous areas according to claim 3, characterized in that, Step S5 specifically includes: Define future projection scenarios, including climate change scenarios and land use change scenarios; Prepare driving data for future scenarios. For climate change scenarios, input future climate prediction data. For land use change scenarios, input future land use maps. Update vegetation cover and management factors and soil and water conservation measures factors in the SWAT model in sync based on future land use types. Run the improved SWAT model to simulate the hydrological, sediment, and phosphorus migration processes in the watershed over a specified future time period; For the phosphorus load of each hydrological response unit output by the improved SWAT model, repeat step S4 to generate a predicted distribution map and trend map of future non-point source phosphorus pollution load.

Images