Phosphorus pollution risk assessment method and system based on migration conversion mechanism

Abstract

This invention relates to the fields of ecological environment and digital simulation, and particularly to a method and system for phosphorus pollution risk assessment based on migration and transformation mechanisms. A method for phosphorus pollution risk assessment based on migration and transformation mechanisms includes the following steps: S10: At preset sampling points in a preset watershed, collect various data according to preset parameter types, wherein the parameter types include the parameter types required by the preset mechanism model and the preset output coefficient model; organize the collected data in chronological order of collection time and establish a spatiotemporally matched database; S20: Extract the parameter types from the database, define the cross-media migration relationship of phosphorus based on physical parameters, construct the transformation relationship of corresponding forms by combining chemical and biological parameters, and classify factor types. This solution solves the problem of low assessment accuracy caused by insufficient sensitivity of static assessment methods and dynamic assessment methods based on simplified models in phosphorus pollution risk assessment.

Description

Technical Field

[0001] This invention relates to the fields of ecological environment and digital simulation, and in particular to a method and system for assessing phosphorus pollution risk based on migration and transformation mechanisms. Background Technology

[0002] In the field of ecological environment, phosphorus pollution risk assessment is the core basis for watershed pollution prevention and control, wastewater treatment plant process optimization, and ecological restoration engineering design. Its accuracy is directly related to the effectiveness of eutrophication control and drinking water source safety, and is a prerequisite for achieving source reduction, process interruption, and end-of-pipe treatment. At present, phosphorus pollution risk assessment is mainly achieved through two types of methods: one is static assessment based on fixed-point monitoring data, using methods such as the single-factor index method and the Nemerow index method, which compares the concentration of total phosphorus and speciation with environmental standards to determine the degree of pollution and solve the basic problem of qualitative assessment of the current pollution status; the other is dynamic assessment based on simplified models, such as using the SWAT model to simulate the migration of phosphorus load in the watershed, or using adsorption and desorption models to calculate the phosphorus release potential in the soil, which initially achieves a rough prediction of pollution trends. However, in riverine terrains with highly uneven spatial distribution of precipitation and significant topographic relief, the patchy distribution of precipitation leads to substantial spatial and temporal differences in soil phosphorus leaching intensity. Furthermore, since phosphorus is mostly particulate, its migration and leaching processes are more significantly influenced by topography and underlying surface conditions (compared to soluble nitrogen). Dramatic changes in surface water velocity caused by topographic relief (such as a sudden increase in velocity on steep slopes and siltation on gentle slopes) further directly alter phosphorus migration pathways and transformation efficiency. In the above process, static evaluation methods can only reflect the pollution status of discrete nodes and are unable to capture this dynamic spatial and temporal heterogeneity. While dynamic evaluation based on simplified models has simulation capabilities, it often treats migration and transformation as independent processes, and the model parameters rely on regional average empirical values, making it difficult to adapt to the differences in hydrodynamics caused by topography and the load fluctuations induced by precipitation. For example, when phosphorus migrates rapidly with high-speed runoff in steep slope areas, the model does not correct the conversion rate parameters in time. In gentle slope siltation areas, phosphorus enrichment and conversion intensify, but the migration flux is not adjusted synchronously. Ultimately, the analysis results are not sensitive enough to pollution hotspots and sudden load changes, the evaluation results deviate significantly from the actual pollution evolution, and the accuracy is difficult to meet the needs of engineering decision-making.

[0003] The migration and transformation mechanism in environmental engineering digital simulation refers to a quantitative model that couples the spatial movement and morphological change processes of matter. Its implementation requires input of physical parameters (such as terrain slope and water flow velocity), chemical parameters (such as phosphorus speciation and pH), and biological parameters (such as the abundance of functional microorganisms). The migration module calculates cross-media flux and pathways, while the transformation module quantifies the morphological transformation rate. A coupling relationship is then established between the two (e.g., using migration flux to correct the initial transformation concentration and using the transformation rate to adjust the migration coefficient), ultimately outputting spatiotemporal data on the dynamic evolution of pollutants. In special environments with uneven precipitation and undulating terrain, the migration and transformation mechanism can accurately adapt to the heterogeneity of hydrodynamics and load inputs, coupling the linkage effect of migration and transformation in real time. This significantly improves the sensitivity of capturing pollution dynamics and fundamentally solves the problem of insufficient accuracy in traditional methods. Summary of the Invention

[0004] This invention provides a method and system for assessing phosphorus pollution risk based on migration and transformation mechanisms, in order to solve the problem that the analysis results of static assessment methods and dynamic assessment methods based on simplified models are not sensitive enough, resulting in low assessment accuracy.

[0005] To solve the above-mentioned technical problems, this application provides the following technical solution: A phosphorus pollution risk assessment method based on migration and transformation mechanisms includes the following steps: S10: At preset sampling points in a preset watershed, collect various data according to preset parameter types, wherein the parameter types include various parameter types required by the preset mechanism model and the preset output coefficient model; organize the collected data in the order of collection time and establish a spatiotemporal matching database; S20: Extract parameter types from the database, define the cross-media migration relationship of phosphorus based on physical parameters, construct the corresponding form transformation relationship by combining chemical and biological parameters, and classify factor types; construct different migration models based on the migration and transformation relationship between different factor types, combined with the preset mechanism model and preset output coefficient model; input various data from the database into different migration models, calculate the key coefficients of each migration model through error inversion, and complete the model optimization; S30: Based on key coefficients, stored environmental quality standards, and comparative calculation results of different migration models, set phosphorus pollution output load standards and screen phosphorus-containing factor types as risk factor types; set prediction features, which include temporal and spatial features; input the currently collected data into the migration model, output phosphorus pollution data that meets the temporal features, and combine the output load standards and spatial features to process the phosphorus pollution data into phosphorus pollution load calculation features and risk factor numerical distribution features, and use the processing results as prediction results; or predict the time that meets the spatial features based on the spatial features as the compliance time, and use the spatial features of the phosphorus pollution data at the compliance time as the prediction result; S40: Integrate the currently collected data, prediction features, and prediction results to generate risk assessment results.

[0006] Furthermore, a phosphorus pollution risk assessment system based on migration and transformation mechanism, employing a phosphorus pollution risk assessment method based on migration and transformation mechanism. The basic principles and beneficial effects of this scheme are as follows: First, this scheme collects all the parameters required for the preset mechanism model and preset output coefficient model and establishes a spatiotemporal matching database to address the lack of spatiotemporal adaptability of the parameters. Then, based on the parameter migration and transformation relationships, it couples the preset mechanism model and preset output coefficient model to construct a migration model. The preset mechanism model is used to characterize the physical migration and chemical and biological transformation processes of phosphorus between interfaces in soil, water, and sediment. It describes microscopic mechanisms such as adsorption-desorption, mineralization, and sedimentation through rate equations, providing a spatially refined dynamic concentration field. The preset output coefficient model rapidly calculates the total watershed load based on land use, fertilizer application, and topographic precipitation correction coefficients, achieving empirical statistics on macroscopic pollution levels. Both methods optimize key coefficients through error inversion, taking into account the advantages of detailed morphological transformation and empirical correction based on topographic precipitation. This allows them to adapt to special topographic watersheds, addressing the shortcomings of traditional models in adapting to complex terrains in plateau and mountainous areas. Finally, they combine temporal and spatial characteristics to conduct predictions, transforming pollution data into load accounting characteristics and risk factor distributions. This enhances the ability to capture pollution dynamics in areas with uneven precipitation and undulating terrain, overcoming the core problems of insufficient sensitivity and accuracy in traditional methods.

[0007] Therefore, this solution addresses the problem of low accuracy in phosphorus pollution risk assessment caused by insufficient sensitivity of analytical results in static assessment methods and dynamic assessment methods based on simplified models.

[0008] Furthermore, in step S10, the water system orientation and soil type distribution characteristics of the preset watershed are obtained, and preset sampling points are set according to the obtained content; the preset sampling points include soil sampling points for different land use types and different planting types, as well as water sampling points at key nodes of the river channel in the catchment area; in step S20, the parameters of the chemical type include phosphorus speciation data and the content of low molecular weight organic acids and humic acid; the database contains soil chemical composition, mineral composition and surface morphology data obtained by X-ray photoelectron spectroscopy, X-ray diffraction, Fourier transform infrared spectroscopy and scanning electron microscopy.

[0009] This scheme deploys sampling points according to water system and soil type, covering different land uses and key river sections, and simultaneously records phosphorus speciation, organic acid and micromineral data, so that the database directly reflects the spatial relationship between source and sink, avoiding the disconnect between sampling and landform; XPS diffraction and electron microscopy results are used to analyze the microscopic mechanism of phosphorus adsorption in soil, which can identify potential release hotspots in advance, make up for the shortcomings of traditional total phosphorus assessment that ignore speciation differences, and improve the accuracy of risk factor identification and spatial prediction.

[0010] Furthermore, in step S20, the extracted physical, chemical, and biological parameters are paired at the same location and time stamp to form a cross-media parameter group. The cross-media migration relationship of phosphorus is defined by the parameter group, and the phosphorus speciation data in the chemical parameters is used as the core driving variable of the transformation relationship. Combined with the abundance of phosphorus cycling microorganisms in the biological parameters, a phosphorus speciation conversion rate equation is constructed. Subsequently, the cross-media parameter group is simultaneously input into the preset mechanism model and the preset output coefficient model. The spatial concentration field output by the preset mechanism model and the regional total load output by the preset output coefficient model are compared grid by grid on the same spatiotemporal grid to generate a dual-model residual matrix. The adjoint sensitivity algorithm is used to perform error inversion on the dual-model residual matrix. The migration rate coefficient, conversion efficiency coefficient of the preset mechanism model, and the correction coefficient of the preset output coefficient model are updated synchronously at one time to realize the mutual constraint between microscopic processes and macroscopic accounting.

[0011] This scheme pairs physical, chemical, and biological parameters at the same spatiotemporal coordinates to form a cross-media parameter set. With phosphorus speciation as the core driver and microbial abundance as the conversion rate coefficient, it simultaneously inputs pre-set mechanistic models and pre-set output coefficient models. This allows for grid-by-grid comparison between the microscopic spatial concentration field and the macroscopic regional total load, generating a residual matrix. Using an adjoint sensitivity algorithm, this matrix is ​​inverted in one operation, simultaneously updating the migration rate coefficient, conversion efficiency coefficient, and correction coefficient, achieving mutual constraints between the two models. This approach eliminates the spatiotemporal misalignment and model contradictions caused by traditional step-by-step inversion. In plateau and mountainous watersheds, a single calculation can simultaneously reveal steep slope scour intensity and potential release on gentle slopes, uncovering hidden mechanisms such as rapid organophosphorus mineralization triggered by microbial proliferation, significantly improving the reliability of subsequent peak and risk window predictions.

[0012] Furthermore, in step S30, the risk factor type screening is associated with the factor types in S20, with soluble exchangeable phosphorus, aluminum-bound phosphorus, and iron-bound phosphorus as the core chemical risk factors, phosphorus cycling microbial abundance as the core biological risk factor, and soil-water migration flux as the core physical risk factor. The output load standard is set according to the preset watershed ecological function zoning, which includes drinking water source protection areas and agricultural irrigation areas. Each level of standard meets the deviation rate threshold of the preset mechanism model and the preset output coefficient model and is dynamically adjusted according to soil type and water flow velocity. When processing the prediction results, the phosphorus pollution load accounting characteristics are weighted and coupled with the corresponding zoning load standard deviation and the proportion of high-value areas in the risk factor numerical distribution characteristics to form a comprehensive risk index.

[0013] This scheme triggers early warnings by multiplying the load index with the spatial risk level map, achieving coordinated quantification of the degree of load exceeding standards and the spatial concentration of risk. It uses a spatiotemporal overlay algorithm to generate a dynamic risk map, which can intuitively mark the diffusion path and the nodes that meet the standards, avoiding the concealment of spatial transmission risks by traditional single concentration indicators. The early warning information is simultaneously pointed to the most sensitive parameter type, providing a clear entry point for subsequent precise regulation and control, significantly improving the identification ability and response speed of scenarios with high risk but no load exceeding standards in complex watersheds such as plateau and mountainous areas.

[0014] Furthermore, in step S30, the temporal characteristics include seasonal cycles, interannual variations, and extreme precipitation event windows, while the spatial characteristics include the upstream catchment area, the midstream main stream area, the downstream estuary area, and key hydrological nodes such as reservoirs and dams. When predicting based on temporal characteristics, each migration model outputs phosphorus pollution load fluctuation curves at different time scales. The periods exceeding the standard are determined by combining the output load standards, and the corresponding risk factor numerical distributions for those periods are superimposed. When predicting based on spatial characteristics, the load thresholds for each sub-unit are first defined. After each migration model outputs the time series that meet the thresholds, the spatial correlation of risk factors for the sub-units at the time of compliance is extracted as a supplementary dimension to the prediction results.

[0015] By allowing temporal features to cover seasonal cycles, interannual variations, and extreme precipitation event windows, and refining spatial features to upstream catchment areas, midstream main stream areas, downstream estuaries, and key nodes of reservoirs and dams, each migration model can output continuous fluctuation curves during temporal prediction. Combined with load standards, the periods exceeding the standard are locked in real time, and the corresponding risk factor distribution is superimposed. During spatial prediction, a load threshold is first set for each sub-unit, and then the spatial correlation of risk factors at the time of compliance is extracted. Thus, the dynamic changes of the periods exceeding the standard and the spatial coupling relationship of key nodes are simultaneously incorporated into the prediction results, forming a risk map with two dimensions of time and space. This significantly improves the accuracy of identifying phosphorus pollution diffusion paths and recovery cycles under conditions of uneven precipitation and topographic relief in plateau and mountainous areas, avoiding misjudgments and omissions caused by the fragmentation of time and space scales in traditional methods.

[0016] Furthermore, in step S30, the phosphorus pollution load accounting characteristics are transformed into a load index based on the output load standard, and the risk factor numerical distribution characteristics are transformed into a three-level spatial risk level map of high, medium and low. The two types of characteristics are used to generate a dynamic risk map through a spatiotemporal overlay algorithm. The map marks the risk diffusion path that meets the time characteristics and the time node that meets the spatial characteristics. A warning threshold is set. When the product of the load index and the area ratio of the high-risk level area exceeds the threshold, a risk warning is automatically triggered. The warning information includes the parameter type with the highest sensitivity among the key coefficients.

[0017] In step S30, a dynamic risk map is generated by spatiotemporally overlaying the load index with the three-level risk level map, and the diffusion path and compliance nodes are marked, realizing the visual fusion of the period of exceeding the standard and the high-risk space. After setting the warning threshold, the automatic warning is triggered by the product of the load index and the area ratio of the high-risk area, and it points to the parameter type with the highest sensitivity. This allows the constructed dual-model optimization system, key coefficient update mechanism and spatiotemporal prediction logic to be presented intuitively, which significantly improves the speed of identifying scenarios where the load has not exceeded the standard but the risk is concentrated and hidden. This ensures that in plateau and mountainous watersheds with uneven precipitation and undulating terrain, the risk diffusion path and recovery time node can be locked in advance and accurately traced to the dominant parameter, avoiding the response lag caused by the fragmentation of indicators in traditional evaluation.

[0018] Further, in step S10, remote sensing data of the remote sensing area is collected using remote sensing equipment. Vegetation data is extracted from the remote sensing images and processed into vegetation indices. The normalized vegetation indices are used to invert the pre-defined dominant plant biomass and phosphorus enrichment. Hyperspectral remote sensing is used to identify phosphorus levels in plant leaves. Simultaneously, infrared remote sensing is used to monitor benthic snail density, waterbird activity, and fecal phosphorus input. The obtained plant phosphorus enrichment and animal phosphorus input are incorporated into a spatiotemporal matching database according to spatial coordinates and timestamps. This database is then used to verify the complementary relationship between biological and chemical processes by calibrating the local migration rate coefficient and conversion efficiency coefficient obtained from the pre-defined mechanism model. In step S20, the soil plant retention and release relationship is defined by plant phosphorus enrichment. The biological sedimentation and transformation relationship in water bodies is constructed by animal population density and fecal phosphorus input. The pre-defined mechanism model adds a plant phosphorus retention module and sets a phosphorus saturation absorption threshold. The pre-defined output coefficient model incorporates animal phosphorus input coefficients. Error inversion simultaneously optimizes the plant phosphorus enrichment coefficient and animal phosphorus deposition coefficient, improving the accuracy of the biologically mediated phosphorus cycle characterization and compensating for the fragmented description of biological processes in traditional models.

[0019] This scheme utilizes remote sensing technology to collect data related to plants (dominant plant biomass, phosphorus accumulation, and leaf phosphorus levels) and animals (density of benthic snails, waterbird activity, and phosphorus input in feces), constructing a complementary verification system for biological and chemical processes. On one hand, it correlates biological data with migration rate coefficients and conversion efficiency coefficients in pre-defined mechanistic models, overcoming the limitations of traditional models that only focus on physicochemical processes. On the other hand, by adding a plant phosphorus interception module, setting a phosphorus saturation absorption threshold, and optimizing plant phosphorus accumulation coefficients and animal phosphorus deposition coefficients, it accurately characterizes the dynamic relationships between soil and plant interception and release, and between water bodies and biological deposition and transformation. This solves the problem of fragmented descriptions of biological participation in the phosphorus cycle in traditional models, making the phosphorus cycle simulation more closely aligned with the synergistic nature of biological and chemical processes in the ecosystem.

[0020] Furthermore, in step S10, special collection is added at the preset sampling points. For human activities, the intensity of farmland fertilization, the frequency of phosphate mining, and the amount of domestic sewage discharge are collected and corresponding interference intensities are set. For microorganisms, key soil layer samples for phosphorus loss are collected. The abundance of phosphate-solubilizing bacteria and polyphosphate-accumulating bacteria is determined by high-throughput sequencing and the phosphatase activity is detected by an enzyme activity analyzer. The interference intensity corresponding to human activities, microbial abundance, and enzyme activity data are entered into a spatiotemporal matching database. In step S20, the human activity data is associated with a preset output coefficient model to correct the phosphorus loss coefficient under different interferences. The microbial data is input into a preset mechanism model and an organophosphate mineralization rate equation is constructed using microbial abundance and enzyme activity. A new human activity phosphorus release coefficient and microbial mineralization efficiency coefficient are added. Through error inversion and synchronous calibration with plant phosphorus enrichment coefficient and animal phosphorus deposition coefficient, a dual-dimensional correction of source strength identification and transformation simulation is achieved, enhancing the ability to characterize phosphorus migration and transformation in phosphate mining areas and farmland composite areas.

[0021] Building upon biological data, this approach incorporates data collection on anthropogenic activities (e.g., farmland fertilization, phosphate mining, domestic sewage) and microorganisms (e.g., abundance of phosphate-solubilizing bacteria and polyphosphate-accumulating bacteria, phosphatase activity), achieving multi-element coverage combining biological, anthropogenic, and microbiological data. By correlating anthropogenic activity data with a pre-defined output coefficient model to correct the phosphorus loss coefficient, and inputting microbial data into a pre-defined mechanism model to construct an organophosphate mineralization rate equation, this approach adds anthropogenic phosphorus release coefficients and microbial mineralization efficiency coefficients, calibrating them synchronously with biological correlation coefficients. This forms a two-dimensional correction mechanism of source strength identification (e.g., anthropogenic interference) and transformation simulation (microbial action). This method specifically addresses the shortcomings of traditional models in characterizing the implicit effects of anthropogenic interference and the microbial-mediated phosphorus transformation process, significantly improving the ability to characterize the complex phosphorus migration and transformation processes in phosphate mining areas and farmland mixed pollution zones. Further, in step S10, remote sensing plant phosphorus enrichment data, animal phosphorus input data, and data on the intensity of human activity disturbance at sampling points, as well as microbial abundance and enzyme activity data, are paired according to the same location and time stamp to form a multi-dimensional parameter set of biological and anthropogenic microorganisms. This set is then incorporated into a spatiotemporal matching database and cross-validated with soil micro-characterization data. In step S20, the multi-dimensional parameter set is simultaneously input into a preset mechanism model and a preset output coefficient model. The preset mechanism model optimizes the adsorption-desorption-mineralization coupling parameters based on plant microbial data, while the preset output coefficient model corrects the regional phosphorus load calculation formula by combining animal and anthropogenic data. When generating the dual-model residual matrix, the residuals of plant-microbial synergistic effects and the residuals of the interaction between anthropogenic mining and animal phosphorus deposition are incorporated. The adjoint sensitivity algorithm is used to synchronously update the plant-microbial coupled phosphorus conversion coefficient and the anthropogenic-animal interaction phosphorus migration coefficient at one time.

[0022] This approach integrates biological, anthropogenic, and microbial data into a multi-dimensional parameter set based on the same location and time stamp. This set is then cross-validated with soil micro-characterization data to construct a data source linking macroscopic processes and microscopic mechanisms. The parameter set is then input into a dual model (a preset output coefficient model and a preset mechanism model) to optimize adsorption-desorption-mineralization coupling parameters, revise the regional phosphorus load calculation formula, incorporate synergistic and interaction residuals, and employ an adjoint sensitivity algorithm to synchronously update coupling and interaction coefficients. This method overcomes the limitations of isolated element characterization, achieving comprehensive synergistic constraints involving biological processes, anthropogenic disturbances, and microbial mediation. The coupled model is upgraded from single-process simulation to system-level coupled simulation, reducing the uncertainty of independent multi-element calculations and strengthening mechanistic support through micro-data validation. Ultimately, this significantly improves the accuracy and reliability of phosphorus pollution risk assessment in complex watershed scenarios. Attached Figure Description

[0023] Figure 1 This is a flowchart of a phosphorus pollution risk assessment method based on migration and transformation mechanism in Example 1. Detailed Implementation

[0024] The following will describe the concept and technical effects of the present invention clearly and completely with reference to embodiments, so as to fully understand the purpose, features and effects of the present invention. Obviously, the described embodiments are only some embodiments of the present invention, not all embodiments. Other embodiments obtained by those skilled in the art based on the embodiments of the present invention without creative effort are all within the scope of protection of the present invention. Example 1

[0025] like Figure 1 As shown, a phosphorus pollution risk assessment method based on migration and transformation mechanisms includes the following steps: S10: Collect various data at preset sampling points in a preset watershed (set by the administrator) according to preset parameter types. The parameter types include various parameter types required by the preset mechanism model and the preset output coefficient model. Organize the collected data in the order of collection time and establish a spatiotemporal matching database. S20: Extract parameter types from the database, define the cross-media migration relationship of phosphorus based on physical parameters, construct the corresponding form transformation relationship by combining chemical and biological parameters, and classify factor types; construct different migration models based on the migration and transformation relationship between different factor types, combined with the preset mechanism model and preset output coefficient model; input various data from the database into different migration models, calculate the key coefficients of each migration model through error inversion, and complete the model optimization; S30: Based on key coefficients, stored environmental quality standards, and comparative calculation results of different migration models, set phosphorus pollution output load standards and screen phosphorus-containing factor types as risk factor types; set prediction features, which include temporal and spatial features; input the currently collected data into the migration model, output phosphorus pollution data that meets the temporal features, and combine the output load standards and spatial features to process the phosphorus pollution data into phosphorus pollution load calculation features and risk factor numerical distribution features, and use the processing results as prediction results; or predict the time that meets the spatial features based on the spatial features as the compliance time, and use the spatial features of the phosphorus pollution data at the compliance time as the prediction result; S40: Integrate the currently collected data, prediction features, and prediction results to generate risk assessment results.

[0026] In step S10, the water system orientation and soil type distribution characteristics of the preset watershed are obtained, and preset sampling points are set according to the obtained information. The preset sampling points include soil sampling points for different land use types and different planting types, as well as water sampling points at key nodes of the river channel in the catchment area. In step S20, the chemical type parameters include phosphorus speciation data and the content of low-molecular-weight organic acids and humic acid; the database contains soil chemical composition, mineral composition, and surface morphology data obtained by X-ray photoelectron spectroscopy, X-ray diffraction, Fourier transform infrared spectroscopy, and scanning electron microscopy.

[0027] Specifically, in step S10, the water system orientation and soil type distribution characteristics of the preset watershed are obtained, and preset sampling points are set based on the obtained information. The preset sampling points include soil sampling points for different land use types and different planting types, as well as water sampling points at key nodes in the catchment area. The preset parameter types specifically include physical, chemical, and biological parameters required for the preset mechanism model. Physical parameters include soil porosity, water flow velocity, and topographic slope; chemical parameters include soluble exchangeable phosphorus, aluminum-bound phosphorus, iron-bound phosphorus, closed-cell phosphorus, calcium-bound phosphorus, residual phosphorus, and total soil phosphorus. The data includes the content of low-molecular-weight organic acids and humic acids, and biological parameters such as the abundance of phosphorus cycling microorganisms. It also includes pollution source characteristic parameters required by the preset output coefficient model, including the area of ​​different land use types, the number of livestock and poultry raised, the population size, the amount of fertilizer applied, and the depth of flooding. During the collection process, the spatial coordinates and sampling time of each sampling point are recorded simultaneously. When organizing the data, it is associated with the corresponding spatiotemporal dimension of the spatiotemporally matched database. The database must contain soil chemical composition, mineral composition, and surface morphology data characterized by X-ray photoelectron spectroscopy, X-ray diffraction, Fourier transform infrared spectroscopy, and scanning electron microscopy.

[0028] During data collection, the location of the water system and the soil type information corresponding to each sampling point need to be recorded simultaneously, and the data should be linked to the corresponding spatial coordinates in the spatiotemporal matching database when processing the data.

[0029] The preset output coefficient model is based on the output coefficient method, incorporating topographic slope correction coefficients extracted from digital elevation models and runoff coefficients calculated from annual average runoff depth and annual precipitation. Its core parameter, the output coefficient, is derived from measured data on different land use types, livestock numbers, and population sizes in a spatiotemporal matching database. Calibration is completed by combining phosphorus loss coefficients and seepage coefficients obtained from watershed phosphorus migration and transformation mechanism studies. The preset mechanism model was selected by comparing its suitability with non-point source pollution simulation models such as SWAT and HSPF. Parameter optimization employs a global sensitivity analysis method. Input data from the spatiotemporal matching database, characterized by soil chemical composition, mineral composition, and surface morphology using X-ray photoelectron spectroscopy, X-ray diffraction, Fourier transform infrared spectroscopy, and scanning electron microscopy, as well as data on phosphorus variation along the course of water bodies, suspended particulate matter, and riverbed sediments, are used to complete the calibration of the localized phosphorus transport coefficient. The localized phosphorus transport coefficient includes the organic phosphorus enrichment ratio in sediment transport and the available phosphorus exchange rate constant in sediments and overlying water. During the simulation process, the two models (output coefficient model and preset mechanism model) share basic meteorological, hydrological and pollution source data in the spatiotemporal matching database. Through error inversion (calculating the residual (error) between the simulated value of the output coefficient model or preset mechanism model and the measured data to clarify the magnitude and direction of the deviation; based on the residual law, targeted adjustment of sensitive parameters (such as migration rate coefficient and mineralization rate constant; repeated iteration until the error meets the threshold to ensure that the model fits the actual phosphorus migration and conversion law), the correction coefficient of the preset output coefficient model and the migration rate coefficient and conversion efficiency coefficient of the preset mechanism model are calculated respectively. The two types of coefficients together serve as the key coefficients corresponding to the migration model.

[0030] In step S20, the extracted physical, chemical, and biological parameters are paired at the same location and time stamp to form a cross-media parameter group. The cross-media migration relationship of phosphorus is defined by the parameter group, and the phosphorus speciation data in the chemical parameters is used as the core driving variable of the transformation relationship. The phosphorus speciation conversion rate equation is constructed by combining the phosphorus cycling microbial abundance in the biological parameters. Then, the cross-media parameter group is simultaneously input into the preset mechanism model and the preset output coefficient model. The spatial concentration field output by the preset mechanism model and the regional total load output by the preset output coefficient model are compared grid by grid on the same spatiotemporal grid to generate the dual-model (output coefficient model and preset mechanism model) residual matrix. The adjoint sensitivity algorithm is used to perform error inversion on the dual-model residual matrix. The migration rate coefficient and conversion efficiency coefficient of the preset mechanism model and the correction coefficient of the preset output coefficient model are updated synchronously at one time to realize the mutual constraint between micro-process and macro-accounting.

[0031] In step S30, the risk factor type screening is associated with the factor types in S20. Soluble exchangeable phosphorus, aluminum-bound phosphorus, and iron-bound phosphorus are the core chemical risk factors, phosphorus cycling microbial abundance is the core biological risk factor, and soil-water migration flux is the core physical risk factor. The output load standard is set according to the preset watershed ecological function zoning, which includes drinking water source protection areas and agricultural irrigation areas. Each level of standard meets the deviation rate threshold of the preset mechanism model and the preset output coefficient model and is dynamically adjusted according to soil type and water flow velocity. When processing the prediction results, the phosphorus pollution load accounting characteristics are weighted and coupled with the proportion of high-value areas (areas exceeding the high-value threshold set by the administrator) in the corresponding zoning load standard deviation and risk factor value distribution characteristics to form a comprehensive risk index.

[0032] In step S30, based on key coefficients, stored environmental quality standards, and comparative calculation results of different migration models, compliance features (including phosphorus migration flux, administrator-preset flow velocity adaptation threshold, administrator-preset phosphorus form concentration threshold, phosphorus form proportion range, phosphorus cycling microbial abundance range, zonal load benchmark, and administrator-preset early warning association threshold) are extracted from the characteristic factor technical map of the compliance watershed input by the administrator. Phosphorus pollution output load standards are set in conjunction with these compliance features, and phosphorus-containing factor types are selected as risk factor types. Prediction features are set, including temporal and spatial features. The currently collected data is input into the migration model, which outputs phosphorus pollution data that meets the temporal features. Combining the set phosphorus pollution output load standards, spatial features, and extracted compliance features, the phosphorus pollution data is processed into phosphorus pollution load calculation features and risk factor numerical distribution features. The processing results are used as prediction results.

[0033] The phosphorus form proportion range is based on chemical risk factors such as soluble exchangeable phosphorus, aluminum-bound phosphorus, and iron-bound phosphorus. It is determined by combining the characteristics of the target watershed with a reasonable range for the proportion of each phosphorus form, which must reduce the risk of phosphorus release (e.g., the proportion of occluded phosphorus and calcium-bound phosphorus is higher than the original management benchmark, while the proportion of soluble exchangeable phosphorus set by the administrator is lower than the safety threshold), and must be compatible with the pre-set watershed ecological function zoning requirements. The phosphorus cycling microbial abundance range refers to the abundance range of phosphorus cycling microorganisms (including phosphate-solubilizing bacteria and polyphosphate-accumulating bacteria) that maintain phosphorus conversion balance. The zoning load benchmark is the phosphorus pollution output load compliance benchmark set by the pre-set watershed ecological function zoning classification.

[0034] In step S30, temporal characteristics include seasonal cycles, interannual variations, and windows of extreme precipitation events, while spatial characteristics include the upstream catchment area, the midstream main stream area, the downstream estuary area, and key hydrological nodes such as reservoirs and dams. When predicting based on temporal characteristics, each migration model outputs phosphorus pollution load fluctuation curves at different time scales. The periods exceeding standards are determined by combining the output load standards with the corresponding risk factor numerical distributions. When predicting based on spatial characteristics, load thresholds for each sub-unit are first defined. After each migration model outputs time series that meet the thresholds, the spatial correlation of risk factors at the compliance time points is extracted as a supplementary dimension to the prediction results.

[0035] In step S30, the phosphorus pollution load accounting characteristics are transformed into a load index based on the output load standard, and the risk factor numerical distribution characteristics are transformed into a three-level spatial risk level map of high, medium, and low. The two types of characteristics are used to generate a dynamic risk map through a spatiotemporal overlay algorithm. The map marks the risk diffusion path that meets the time characteristics and the time node that meets the spatial characteristics. An early warning threshold is set. When the product of the load index and the area ratio of the high-risk level area (the load index only quantifies the phosphorus pollution intensity of the preset watershed, and the area ratio of the high-risk level area only reflects the spatial range of pollution, and a single indicator is prone to misjudgment; the product of the two integrates the preset watershed phosphorus pollution intensity and the spatial range of pollution to accurately quantify the overall threat and avoid one-sided assessment) exceeds the threshold (the administrator sets the upper limit based on the empirical data of the compliant watershed, combines it with the environmental quality standard calibration, verifies it through dual models and adjusts it according to ecological zoning; it is used to judge whether the product of the load index and the area ratio of the high-risk area exceeds the critical value, and assists in screening risk factor limits, providing a standard for judging whether the pollution exceeds the carrying capacity and whether regulation is needed) the risk warning is automatically triggered. The warning information includes the parameter type with the highest sensitivity among the key coefficients.

[0036] Specifically, the system acquires basic meteorological data (including spatial variation data of pre-defined watershed patch precipitation input by the administrator, such as monthly precipitation in the southeast region being more than three times that of the northwest region during the rainy season; extreme precipitation window data set by the administrator according to stored standards, such as defined as daily precipitation ≥50mm and lasting ≥6 hours) and basic pollution source data (including annual fertilizer application in farmland; annual mining frequency of phosphate mines, counted by the number of mining operations per year). It categorizes the data by seasonal cycle into rainy season or dry season (based on the precipitation threshold set by the administrator, and binding precipitation intensity correction coefficients to the rainy and dry seasons respectively, for example, 1.3-1.6 for the rainy season and 0.6-0.8 for the dry season). It also categorizes the data by interannual pollution source variation into high, medium, and low load years (based on the pollution source threshold set by the administrator, and calibrating the organophosphate mineralization rate constant, for example, 0.05d⁻¹ for high load years and 0.02d⁻ for low load years). ¹) The extreme precipitation window period is divided into three sub-periods. The adjustment methods for the migration rate coefficient, precipitation intensity, and precipitation intensity weight of the three sub-periods and different sub-periods are all based on the dynamic changes of hydrological, chemical, and biological processes during extreme precipitation in the watershed (assuming that it is set to 12 hours before precipitation in the pre-leaching period, 6-12 hours after precipitation in the strong migration period, and 12-24 hours after precipitation in the delayed sedimentation period, and the migration rate coefficient is adjusted for different sub-periods, such as taking the coefficient of 0.05 cm / s when the flow velocity on the steep slope section in the strong migration period is >1.5 m / s). Quantitative values ​​are assigned to each time sub-feature (such as assigning values ​​of 3, 2, and 1 to the three levels of extreme, heavy rain, and light rain, with the weight increased by 20% in the rainy season and decreased by 50% in the dry season). Finally, the time type and parameter correction values ​​are integrated into a feature vector containing dimensions such as time type (such as 0 for the pre-leaching period, 1 for the strong migration period, and 2 for the delayed sedimentation period) and precipitation level.

[0037] Based on stored technical standards, preset watershed measured data, and dual-model adaptation requirements, the administrator sets slope grading, land use assignment, flow velocity gradient, and soil adsorption coefficient gradient, taking into account the soil mineral composition patterns associated with topography, the differences in pollution source intensity among different land uses, and the dynamic response characteristics of river flow velocity. The administrator also determines the grid cell size based on the requirements for capturing watershed spatial heterogeneity, computational efficiency, and compatibility with external monitoring data. Finally, referencing the measured data of the superposition effect of dual-model coupling and the results of risk sensitivity analysis, the administrator sets correction parameters for high-risk cells, thus adapting the overall model error threshold and pollution targeted control requirements.

[0038] First, acquire the pre-defined watershed DEM (Digital Elevation Model) topographic data (assuming steep slopes >25°, gentle slopes 5°-25°, and flat slopes <5°), land use data (assuming steep slope farmland, gentle slope forest, etc.), river monitoring data (such as flow velocity and sedimentation rate at each section), and soil mineral data (such as high quartz on steep slopes and high iron-aluminum oxide on gentle slopes); divide the watershed into grid cells of a pre-defined size (assuming 1km × 1km), and assign parameters to slope (assuming a spatial parameter of 3 for steep slopes, 2 for gentle slopes, and 1 for flat slopes), land use (… Assuming that the spatial parameters of cultivated land are 3, forest land is 1, and grassland is 2, and that the river flow velocity is quantified (assuming that the spatial parameters >1.5m / s are 3, 0.3-1.5m / s are 2, and <0.3m / s are 1), the spatial parameters of high-risk units such as steep-slope cultivated land and steep-slope river channels are corrected (assuming that the product of the slope and flow velocity is increased by 1.5 times). Combined with the soil adsorption coefficient (assuming 0.1 for steep slopes and 0.8 for gentle slopes), the grid cells and corresponding spatial parameters are integrated into a 1×7-dimensional spatial feature vector to support the grid-by-grid comparison and calculation of the two models.

[0039] In practice, in river basins with extremely uneven spatial distribution of precipitation and large topographic relief, the slope classification extracted by the basin digital elevation model (DEM) (for example, divided into three levels, each corresponding to steep slopes >25°, gentle slopes 5°-25°, and flat slopes <5°) is combined with the patchy precipitation data from precipitation monitoring stations (assuming that the monthly precipitation in the southeastern part of the basin during the rainy season is more than 3 times that in the northwestern part).

[0040] Soil sampling points were set up in steep-slope farmland (such as tobacco planting areas), gentle-slope woodland, and flat-slope river valley grassland. Simultaneously, water sampling points were set up in steep-slope river sections (e.g., flow velocity > 1.5 m / s), gentle-slope river sections (e.g., flow velocity 0.3-1.5 m / s), and estuary siltation areas. When collecting parameters, the instantaneous precipitation intensity, soil moisture content (to reflect phosphorus leaching potential), and river flow velocity were recorded at each point. Linkage data between slope, flow velocity, and phosphorus speciation were also added (assuming that due to rapid runoff inhibiting adsorption, the proportion of soluble exchangeable phosphorus in the soil of steep slopes is twice that of gentle slopes). The database should include mineral composition differences between soils with high quartz content (reducing phosphorus adsorption capacity) on steep slopes and those with high iron and aluminum oxide content (enhancing phosphorus adsorption) on gentle slopes, as analyzed by X-ray diffraction.

[0041] When constructing a migration model, the cross-media parameter set needs to be labeled with precipitation intensity, slope, and flow velocity. For example, for steep slopes, physical parameters (including flow velocity and slope) are paired with chemical parameters (including soluble exchangeable phosphorus) and biological parameters (including phosphorus cycling microbial abundance, which is low on steep slopes due to strong soil erosion) to define the rapid phosphorus migration relationship driven by high-velocity runoff; for gentle slopes, physical parameters (e.g., sedimentation rate) are associated with chemical parameters (e.g., iron-bound phosphorus) to construct the phosphorus transformation relationship under sedimentation conditions.

[0042] The preset mechanistic model focuses on simulating the non-equilibrium migration of phosphorus with high-velocity runoff on steep slopes (e.g., adjusting the adsorption rate by inputting quartz mineral content data) and the microbial transformation of phosphorus in sediments on gentle slopes (e.g., optimizing the desorption constant by inputting iron and aluminum oxide data). The preset output coefficient model sets differentiated loss coefficients for different slope zones (assuming the coefficient for arable land on steep slopes is 5 times that of forest land on gentle slopes). It is assumed, through the dual-model residual matrix, that the deviation between the simulated values ​​of the mechanistic model and the calculated values ​​of the output coefficient model on steep slopes mainly stems from insufficient phosphorus desorption caused by a sudden increase in flow velocity. Based on this, the migration rate coefficient of the preset mechanistic model and the slope correction coefficient of the preset output coefficient model are adjusted simultaneously.

[0043] In forecasting, the temporal characteristics are focused on the extreme precipitation window during the rainy season (assuming a single-day precipitation exceeding 50 mm), while the spatial characteristics are refined to steep slope confluence areas, gentle slope sedimentation zones, and estuary diffusion areas. Based on the temporal characteristics, the migration model outputs a curve showing a sudden increase in soluble phosphorus load on steep slopes within a period after precipitation (assumed to be within 24 hours) (assumed to exceed the storage standard for agricultural irrigation areas by 30%). Overlaying this curve shows that the high-risk area during this period is concentrated between steep slope farmland and a downstream gentle slope river channel at a distance (assumed to be 1 km).

[0044] When predicting based on spatial characteristics, after defining the load threshold for steep slope sections, it was found that the compliance time point was extracted and found that gentle slope sections require 7 days longer to meet the standard due to phosphorus deposition than steep slope sections, and the spatial correlation of risk factors between the two reached 0.8 (indicating that steep slope output is the main source of risk for gentle slopes). The dynamic risk map marks the risk paths of precipitation, steep slope erosion, gentle slope deposition, and estuary diffusion. When the product of the rainy season load index and the proportion of the high-risk area on the steep slope exceeds the warning threshold (specifically set by the administrator), an early warning is triggered, and the velocity sensitivity coefficient is indicated as a key control parameter (for example, for every 0.1 m / s increase in velocity on a steep slope section, phosphorus migration flux increases by 15%).

[0045] Example 2 The only difference between this embodiment and Embodiment 1 is that, in step S10, sampling points are set up in the catchment area based on the watershed's river system orientation, soil type distribution (e.g., mainly red and yellow soil), and pollution source distribution (including industry, urban life, tourism, agricultural non-point source pollution, and free-range livestock farming). Soil sampling points cover cultivated land (e.g., tobacco and wheat planting areas), grassland, forest land, and the soil surrounding phosphate mining enterprises, phosphate chemical production enterprises (e.g., phosphate fertilizer production, yellow phosphorus preparation), and phosphogypsum storage enterprises. Simultaneously, gradient diffusion film technology (DGT) is used to construct a soil available phosphorus pool, and the migration kinetics parameters of the soil available phosphorus interface are collected. Water sampling points include key river nodes (e.g., tributary confluences, downstream of the discharge outlets of phosphate mining enterprises, phosphate chemical production enterprises, and phosphogypsum storage enterprises). Simultaneously, water, suspended particulate matter, and riverbed sediment samples were collected from reservoirs and dams. Soil organic phosphorus active component data were collected using a modified Bowman-Cole organic phosphorus classification method. A soil column culture experiment was designed to simulate soil phosphorus leaching patterns and seepage coefficients under different rainfall intensities (assuming mild rainfall of <10 mm / d, moderate rainfall of 10-25 mm / d, and severe rainfall of >25 mm / d) and fertilization amounts (e.g., conventional fertilizer and reduced-volume fertilizer; conventional fertilizer refers to the amount and ratio of fertilizer commonly used by local farmers or planting entities for specific crops in agricultural production in the study area; reduced-volume fertilizer refers to a fertilization gradient formed by reducing the amount of phosphorus applied on the basis of "conventional fertilizer," with the amount applied being lower than the conventional fertilization level, and the specific reduction ratio determined in combination with the phosphorus pollution control needs and crop growth needs of the study area).

[0046] In addition to physical parameters (including soil porosity, water flow velocity, and topographic slope, such as three slope classifications extracted from DEM: steep slope >25°, gentle slope 5-25°, and flat slope <5°), chemical parameters (including phosphorus forms such as S / LP, Al-P, and Fe-P, and low-molecular-weight organic acids), and biological parameters (including phosphorus cycling microbial abundance), the database adds migration factor data such as soil erosion modulus, sub-basin shape coefficient, and distance to rivers, as well as pollution source data such as industrial pollution emissions and tourism wastewater generation. The database is linked to the spatial coordinates of sampling points, sampling time, and experimental data, including soil mineral data characterized by XPS and XRD, available phosphorus data determined by DGT, and seepage coefficient data from soil column experiments, achieving spatiotemporal matching.

[0047] In step S20, phosphorus migration pathways are divided into soil and water migration factors (e.g., related to topographic slope and soil erosion modulus), sediment and water exchange factors (e.g., related to the phosphorus exchange rate constant between sediment and overlying water, specifically measured by indoor culture experiments), and pollution source output factors (e.g., related to the intensity of industrial and livestock breeding pollution sources).

[0048] The output coefficient model is improved to obtain the preset output coefficient model: In view of the characteristics of the plateau monsoon climate, uneven spatial distribution of precipitation, and large topographic relief of the watershed, the topographic slope correction coefficient (Kslope) and precipitation intensity correction coefficient (Kpre) are added to the traditional output coefficient model. The adjustment formula is shown in the following formula (1): (1), Where L represents the total phosphorus loss in the watershed (unit: kg / a). This is a terrain slope correction factor (assuming a value of 1.2-1.5 for steep slopes, 0.8-1.1 for gentle slopes, and 0.5-0.7 for level slopes, with the specific value determined based on the DEM slope classification). The precipitation intensity correction factor is assumed to be 1.3-1.6 for heavy rainfall, 0.9-1.2 for moderate rainfall, and 0.6-0.8 for light rainfall, with the specific value determined based on watershed precipitation monitoring data. The phosphorus export coefficients (unit: kg / (hm²・a)) for different land use types were obtained by constructing a localized phosphorus transport coefficient set, based on soil column culture experiment data, watershed field monitoring load data, and combined with stored regional similar coefficients, and then corrected through dual-model comparison and verification. Land use area (unit: hm²); The phosphorus input in the land use area is expressed in kg / hm². For farmland areas, the data is derived from a combination of field monitoring (weighing of fertilizer application) and farmer questionnaires. For other areas, the phosphorus background value from the soil survey is used as a reference. The phosphorus output coefficient for livestock and poultry farming (unit: kg / (head·a) or kg / (animal·a)) is obtained by calibrating a preset mechanism model and a preset output coefficient model based on actual measured data of sewage discharge from livestock and poultry farming and coefficients of similar farming models obtained by the administrator. The number of livestock and poultry raised (unit: head / animal) is determined by on-site verification of the number of animals in the farms and the breeding statistics of the pre-defined watershed. The phosphorus input for livestock and poultry farming (unit: kg / head or kg / animal) is calculated based on the phosphorus content test data of the feed and the feed consumption during the farming cycle. The output coefficient (unit: kg / (person·a)) is derived by back-calculating phosphorus emission data from the monitoring of domestic sewage in the preset watershed, combined with population statistics yearbook data. The population of the basin (in people) is obtained by using the resident population statistics of the preset basin. The per capita phosphorus input (unit: kg / person) is calculated based on regional residents' dietary structure survey and phosphorus content testing data of daily consumer products. The phosphorus output coefficient of industrial pollution sources (unit: kg / m³) is derived from on-site monitoring data of industrial enterprise discharge outlets and enterprise discharge declaration data, and is obtained after verification and confirmation. Industrial wastewater discharge (unit: m³ / a) is obtained from the acquired enterprise wastewater discharge data; P is the total phosphorus input from rainfall (unit: kg / a), which is obtained by using preset watershed precipitation monitoring data and phosphorus deposition sampling analysis to obtain the deposition rate, combined with the watershed area and runoff coefficient calculated by GIS, and calculated based on an empirical statistical model. The calculation method is non-point source phosphorus pollution load accounting analysis based on an empirical statistical model.

[0049] By using dual correction coefficients, the model accurately characterizes the pre-defined watershed topography (plateau monsoon climate, watershed topography with uneven spatial distribution of precipitation and large topographic relief) and the impact of precipitation on phosphorus loss. At the same time, it incorporates industrial pollution source items to solve the problem of insufficient coverage of complex pollution sources in mountainous areas by traditional models.

[0050] In this embodiment, the AnnAGNPS model, suitable for small watersheds in mountainous areas (different from the SWAT / HSPF model in Example 1), was selected. Localized parameters were input, including soil parameters and hydrological parameters. Among them, the soil parameters were the available phosphorus interface migration kinetic parameters measured by DGT and the organic phosphorus active components measured by the modified Bowman-Cole method; the hydrological parameters were the seepage coefficient simulated by the soil column experiment and the sub-watershed shape coefficient measured in the watershed. The global sensitivity analysis method was used to optimize the AnnAGNPS model parameters (such as the phosphorus exchange rate constant between sediment and overlying water and the phosphorus use efficiency of phytoplankton), and the migration rate coefficient and conversion efficiency coefficient of the model were updated synchronously through error inversion.

[0051] In step S30, the output load standard is set according to a preset phosphorus index classification. It is assumed that, based on the phosphorus index model, the risk levels are divided into a pollution-free maintenance zone (PI < 0.3), a low-pollution identification zone (0.3 ≤ PI < 0.6), a medium-pollution concern zone (0.6 ≤ PI < 0.9), and a key pollution warning zone (PI ≥ 0.9). Each standard is associated with threshold values ​​for factors such as soil erosion modulus and phosphate fertilizer application rate.

[0052] The risk factor type is associated with the factor type in step S20, including source factors (e.g., phosphorus environmental threshold, land use index, phosphate fertilizer application rate) and migration factors (e.g., soil erosion modulus, annual runoff depth, sub-basin shape coefficient, distance to river).

[0053] The risk is quantified using an improved phosphorus index formula, as shown in formula (2) below: (2), PI is the phosphorus pollution risk index (dimensionless). The environmental threshold score for phosphorus is assumed to be 1-5 points, with higher scores for lower thresholds. The land use index grade score is assumed to be 1-5 points, with cultivated land score > grassland score > forest land score. The score is based on the application rate of phosphate fertilizer (assuming a score of 1-5, with higher application rates resulting in higher scores). The source factor weights are assumed to be 0.3, 0.4, and 0.3 respectively, and the specific values ​​will be determined through expert consultation. The score is the annual runoff depth rating (assuming a range of 1-5 points, with higher scores for greater runoff depth). The score is the soil erosion modulus level (assuming a score of 1-5, with a higher score for a larger modulus). This is the score for the subsurface drainage level (assuming a score of 1-5, with higher scores for stronger drainage). The distance to the river is assigned a score (assuming a score of 1-5, with higher scores for closer distances). The migration factor weights are assumed to be 0.25, 0.3, 0.2, and 0.25 respectively, and the specific values ​​will be determined through expert consultation.

[0054] Formula (2) strengthens the synergistic effect of source strength (characterizing the intensity or ability of phosphorus pollution sources to release phosphorus into the environment (soil, water bodies)) and migration ability by multiplying the source factor and migration factor, and more accurately characterizes the high-risk characteristics of steep slope farmland with high source strength and short distance to the river in the preset watershed.

[0055] Dynamic risk maps identify pollution sources, migration paths, and risk levels. For example, they depict risk diffusion paths such as upstream livestock and poultry farming areas, leaching from steep slope farmland, and tributary confluence with the main stream. Based on PI values ​​and load calculation results, key source areas are identified, namely upstream tobacco planting areas (assuming PI ≥ 0.9, contributing 40% of the total load) and the areas surrounding "three phosphorus" enterprises (assuming PI ≥ 1.0, contributing 25% of the total load), providing targeted areas for precise prevention and control.

[0056] In step S40, the risk assessment results are integrated by combining the currently collected DGT available phosphorus data and soil column test seepage coefficient data with the phosphorus index prediction results and key source area identification results from step S30 to generate a four-dimensional assessment report that combines data, models, risks, and source areas. This report clarifies the prevention and control recommendations that prioritize the reduction of fertilizer application on upstream farmland and the upgrading of wastewater from "three phosphorus" enterprises, and is adapted to the characteristics of complex pollution in plateau and mountainous areas.

[0057] In practice, taking key polluted river basin A as the target, the implementation will follow the core process as follows: Based on the watershed system, soil type, and pollution sources (including industry, livestock and poultry farming, and phosphorus-producing enterprises), soil sampling points were set up covering cultivated land (e.g., tobacco and wheat planting areas), grassland, woodland, and the vicinity of phosphate mines / phosphate chemical plants / phosphogypsum storage facilities. Soil available phosphorus was measured using gradient diffusion film technology (DGT), and organic phosphorus active components were collected using the modified Bowman-Cole method. Water sampling points included tributary confluences, downstream of phosphorus-producing enterprise discharge outlets, and reservoir dams, simultaneously collecting water samples, suspended particulate matter, and sediment samples. A soil column experiment was designed to simulate, for example, three rainfall intensities (assumed to be <10 mm / d (light), 10-25 mm / d (moderate), and >25 mm / d (severe)) and two fertilizer application rates (conventional fertilizer and reduced fertilizer), to measure soil phosphorus leaching patterns and permeability coefficients.

[0058] The parameters include physical (including soil porosity and slope classification extracted from DEM), chemical (including phosphorus forms such as S / LP and Al-P), and biological (including phosphorus cycling microbial abundance). Newly added migration factors such as soil erosion modulus and sub-basin shape coefficient, as well as industrial pollution source data, are included. The database is linked with sampling point coordinates, time, and experimental data (including soil mineral data characterized by XPS / XRD and DGT available phosphorus data) to achieve spatiotemporal matching.

[0059] Based on migration pathways, factors are categorized as follows: soil and water migration factors (including associated slope and erosion modulus), sediment and water exchange factors (including associated phosphorus exchange rate constants measured in laboratory experiments), and pollution source output factors (including associated industrial or livestock pollution intensity).

[0060] The improved output coefficient model incorporates terrain slope correction coefficients (assuming Kslope values ​​are 1.2-1.5 for steep slopes, 0.8-1.1 for gentle slopes, and 0.5-0.7 for level slopes) and precipitation intensity correction coefficients (assuming Kpre values ​​are 1.3-1.6 for severe, 0.9-1.2 for moderate, and 0.6-0.8 for mild), and includes industrial pollution source terms. The AnnAGNPS model is selected, and localized parameters (including DGT available phosphorus and soil column seepage coefficient) are input. Global sensitivity analysis is used to optimize parameters such as sediment phosphorus exchange rate, and error inversion is used to update the model coefficients.

[0061] The standard is set according to the phosphorus index (assuming no pollution PI<0.3, low pollution 0.3≤PI<0.6, medium pollution 0.6≤PI<0.9, and key warning PI≥0.9), and the erosion modulus and phosphate fertilizer application threshold are associated; the risk factors include source factors (including phosphorus threshold, land use index, and phosphate fertilizer application) and migration factors (including erosion modulus, annual runoff depth, etc.).

[0062] Risk is quantified using an improved phosphorus index formula (e.g., source factor × migration factor), and a dynamic map is used to mark pollution sources and migration paths; key source areas are identified, including upstream tobacco planting areas (assuming PI ≥ 0.9, contributing 40% of the load) and areas surrounding "three phosphorus" enterprises (assuming PI ≥ 1.0, contributing 25%).

[0063] By integrating DGT data, soil column permeability coefficient, phosphorus index results, and key source area information, a four-dimensional report combining data, models, risks, and source areas is generated, which recommends reducing fertilizer application on upstream farmland and upgrading wastewater standards for enterprises producing phosphorus, phosphorus, and phosphorus products, and is adapted to the characteristics of complex pollution in mountainous areas.

[0064] Example 3 The only difference between this embodiment and embodiments 1-2 is that in step S10, remote sensing data of the remote sensing area is collected by remote sensing equipment, vegetation data is extracted from the remote sensing images and processed into vegetation indices, and the pre-defined dominant plant biomass and phosphorus enrichment are inverted by normalizing the vegetation indices; the phosphorus level of plant leaves is identified by hyperspectral remote sensing, and the density of benthic snails, waterbird activity and phosphorus input in feces are monitored simultaneously by infrared remote sensing. The obtained plant phosphorus enrichment and animal phosphorus input are included in the spatiotemporal matching database according to spatial coordinates and timestamps, and the local migration rate coefficient and conversion efficiency coefficient obtained by calibration with the pre-defined mechanism model form a complementary verification of biological and chemical processes.

[0065] In step S20, the relationship between soil plant retention and release is defined by the amount of phosphorus enrichment in plants, and the relationship between biological sedimentation and transformation in water bodies is constructed by animal population density and phosphorus input in feces. The preset mechanism model adds a plant phosphorus retention module and sets a phosphorus saturation absorption threshold. The preset output coefficient model incorporates animal phosphorus input coefficients. Error inversion simultaneously optimizes the plant phosphorus enrichment coefficient and animal phosphorus deposition coefficient, improving the accuracy of bio-mediated phosphorus cycle characterization and making up for the fragmented description of biological processes in traditional models.

[0066] In step S10, special collection is added at the preset sampling points. For human activities, the intensity of farmland fertilization, the frequency of phosphate mining, and the amount of domestic sewage discharge are collected and corresponding interference intensities are set. For microorganisms, key soil layer samples for phosphorus loss are collected. The abundance of phosphate-solubilizing bacteria and polyphosphate-accumulating bacteria is determined by high-throughput sequencing and phosphatase activity is detected by enzyme activity analyzer. The interference intensity corresponding to human activities, microbial abundance, and enzyme activity data are entered into the spatiotemporal matching database.

[0067] In step S20, human activity data is correlated with a preset output coefficient model to correct the phosphorus loss coefficient under different disturbances. Microbial data is input into a preset mechanism model and an organophosphorus mineralization rate equation is constructed using microbial abundance and enzyme activity. A new human activity phosphorus release coefficient and microbial mineralization efficiency coefficient are added. Through error inversion and synchronous calibration with plant phosphorus enrichment coefficient and animal phosphorus deposition coefficient, a dual-dimensional correction of source strength identification and transformation simulation is achieved, enhancing the ability to characterize phosphorus migration and transformation in phosphate mining areas and farmland complex areas.

[0068] In step S10, remote sensing plant phosphorus enrichment data, animal phosphorus input data, and data on human activity disturbance intensity, microbial abundance, and enzyme activity at sampling points are paired according to the same location and time stamp to form a multi-dimensional parameter group of biological and human-made microorganisms. This group is then incorporated into a spatiotemporal matching database and cross-validated with soil micro-characterization data.

[0069] In step S20, the multi-dimensional parameter set is simultaneously input into the preset mechanism model and the preset output coefficient model. The preset mechanism model optimizes the adsorption-desorption-mineralization coupling parameters based on plant microbial data. The preset output coefficient model combines animal and anthropogenic data to correct the regional phosphorus load calculation formula. When generating the dual-model residual matrix, the residual of plant-microbial synergy and the residual of the interaction between anthropogenic mining and animal phosphorus deposition are included. The adjoint sensitivity algorithm is used to update the plant-microbial coupled phosphorus conversion coefficient and the anthropogenic-animal interaction phosphorus migration coefficient at one time.

[0070] In practice, the study area is phosphorus-rich zone A, which is assumed to be approximately 10 km², including 2 km² of phosphate-bearing mining area, 5 km² of farmland, and 3 km² of forest land. The implementation follows the process of data collection, model building, and parameter optimization, as follows: 250m resolution MODIS NDVI data (assuming cloud cover <10%, such as an image from March of a certain year), 10m resolution Sentinel-2 hyperspectral data, and 30m resolution Landsat-8 infrared data were selected, downloaded from the corresponding satellite data platforms, and processed using ENVI 5.6.

[0071] First, the spatial range of vegetation distributed around phosphorus-rich zone A (e.g., Rumex japonica), on barren slopes (Sedum candelilla) in phosphorus-rich zone A, and on the edge of forest land in phosphorus-rich zone A (e.g., Miscanthus sinensis) is delineated through remote sensing interpretation.

[0072] Assuming that three 10m×10m field plots are set up for each type of plant; the mean NDVI of the euphorbia lysate plot is 0.6, the field biomass is 2000kg / hm², and the leaf phosphorus content is 11.5g / kg, the biomass is fitted to be 3000×NDVI+200, and the leaf phosphorus content is 18×(750nm / 550nm)-2, and the phosphorus enrichment of plants in the whole region is inverted.

[0073] Infrared remote sensing is used to identify benthic snail habitats (e.g., shallow water areas with a water temperature of 18℃). Assuming there are three river monitoring sections in phosphorus-rich area A, with a snail density of 80 individuals / m² at each section, the total snail population in the area is calculated. Waterbirds are identified through infrared thermal anomalies; 3000 waterbirds were observed quarterly, and based on an average daily phosphorus input of 0.02g per bird in feces, the quarterly input is calculated to be 1.8kg.

[0074] In phosphorus-rich area A, multiple (e.g., 5) monitoring plots were set up for farmland fertilization. It was assumed that the fertilizer application rate was 300 kg / hm² for conventional fertilizer and 180 kg / hm² for reduced fertilizer, and the results were confirmed through questionnaires and weighing. The phosphate mining area was checked quarterly. It was assumed that mining lasted 15 days in Q1, and the interference intensity was classified as moderate. Three pipeline flow meters were set up in the residential area. It was assumed that the average daily sewage discharge was 80 t, and the interference intensity was classified as moderate.

[0075] Multiple (e.g., 3) sampling points were set up in the phosphate mining area, farmland, and woodland of phosphorus-rich zone A. Soil layers of 0-10cm and 10-20cm were collected, with 3 replicates per layer. High-throughput sequencing was performed using phoD gene primers, assuming an abundance of phosphate-solubilizing bacteria of 10 in the phosphate mining area. 7 copies / g, farmland 8×10 6 Copies / g; Phosphatase activity was measured using the PNPP method, assuming a measured value of 0.8 μmol / (g·h) in phosphate mining areas and 0.4 μmol / (g·h) in forest areas, and the data were entered into a spatiotemporal matching database.

[0076] Assuming a phosphorus-rich zone A is divided into 1km×1km grids, with one core point in each grid, and the following treatments are applied: 0.6% NDVI of *Rumex acetosa*, 300 kg / hm² of farmland fertilizer, and 10% phosphate-solubilizing bacteria. 7 Data such as copies / g were paired with quarterly timestamps and cross-validated with soil organic phosphorus percentages obtained from XPS analysis (assuming 25% organic phosphorus in areas with high phytophosphorus enrichment and 10% in surrounding areas). A phosphorus saturation absorption threshold of 11.5 g / kg was set for *Rumex japonicus* leaves and 7.99 g / kg for *Imperata cylindrica* leaves; absorption calculations ceased when leaf phosphorus content reached this threshold. An organic phosphorus mineralization equation (e.g., dP / dt = k × M × E) was used, assuming k = 0.05 (high abundance of phosphate-solubilizing bacteria) and M = 10 in phosphate mining areas. 7 Given copies / g, E = 0.8 μmol / (g·h), the calculated mineralization rate is 4 × 10⁻⁶. 6 mg / (kg・d).

[0077] Assuming a waterbird population of 0.02 kg / (bird·a) and a snail population of 0.001 kg / (bird·a), and an annual waterbird population of 5000 in phosphorus-rich area A, the annual input is calculated to be 100 kg. The phosphate mine Q1 is mined for 15 days, with a coefficient of 1.2 (mining frequency 50%). The coefficients for conventional fertilizer and reduced-volume fertilizer in farmland are 1.0 and 0.6, respectively. The residual weight for plant-microbe synergy is 0.4 (assuming a simulated mineralization amount of 5 × 10⁻⁶). 6 mg / (kg・d), measured 4×10 6 Residual 1×10 6The residual weights for human-animal interaction were set at 0.6. Using SWAT-CUP software, the initial migration rate coefficient was 0.03 cm / s. After error inversion, the plant-microbe coupling coefficient was updated to 0.04 and the human-animal interaction coefficient was updated to 0.02. The final model had R²=0.88 and NSE=0.82, which met the accuracy requirements.

[0078] The above are merely embodiments of the present invention. The invention is not limited to the fields covered by these embodiments. Commonly known structures and characteristics in the solutions are not described in detail here. Those skilled in the art are aware of all common technical knowledge in the field prior to the application date or priority date, are able to access all existing technologies in that field, and have the ability to apply conventional experimental methods prior to that date. Those skilled in the art can, under the guidance of this application, improve and implement this solution in combination with their own capabilities. Some typical known structures or methods should not be obstacles for those skilled in the art to implement this application. It should be noted that those skilled in the art can make several modifications and improvements without departing from the structure of the present invention. These should also be considered within the scope of protection of the present invention, and will not affect the effectiveness of the implementation of the present invention or the practicality of the patent. The scope of protection claimed in this application should be determined by the content of its claims, and the specific embodiments described in the specification can be used to interpret the content of the claims.

Claims

1. A method for assessing phosphorus pollution risk based on migration and transformation mechanisms, characterized in that, Includes the following steps: S10: At preset sampling points in a preset watershed, collect various data according to preset parameter types, wherein the parameter types include various parameter types required by the preset mechanism model and the preset output coefficient model; organize the collected data in the order of collection time and establish a spatiotemporal matching database; S20: Extract parameter types from the database, define the cross-media migration relationship of phosphorus based on physical parameters, construct the transformation relationship of corresponding forms by combining chemical and biological parameters, and classify factor types; Based on the transfer and transformation relationships between different factor types, different transfer models are constructed by combining a preset mechanism model and a preset output coefficient model; various data from the database are input into different transfer models, and the key coefficients of each transfer model are calculated through error inversion to complete model optimization; S30: Based on key coefficients, stored environmental quality standards, and comparative calculation results of different migration models, set phosphorus pollution output load standards and screen phosphorus-containing factor types as risk factor types; The prediction features are set, which include temporal and spatial features; The currently collected data is input into the migration model, which outputs phosphorus pollution data that meets the time characteristics. Combined with the output load standard and spatial characteristics, the phosphorus pollution data is processed into phosphorus pollution load accounting characteristics and risk factor numerical distribution characteristics. The processing results are used as prediction results; or the time that meets the spatial characteristics is predicted based on the spatial characteristics as the time to meet the standard, and the spatial characteristics of the phosphorus pollution data at the time of meeting the standard are used as the prediction results. S40: Integrate the currently collected data, prediction features, and prediction results to generate risk assessment results.

2. The phosphorus pollution risk assessment method based on migration and transformation mechanism according to claim 1, characterized in that: In step S10, the water system orientation and soil type distribution characteristics of the preset watershed are obtained, and preset sampling points are set according to the obtained information. The preset sampling points include soil sampling points for different land use types and different planting types, as well as water sampling points at key nodes of the river channel in the catchment area. In step S20, the parameters of the chemical type include phosphorus speciation data and the content of low molecular weight organic acids and humic acid; the database contains soil chemical composition, mineral composition and surface morphology data obtained by X-ray photoelectron spectroscopy, X-ray diffraction, Fourier transform infrared spectroscopy and scanning electron microscopy.

3. The phosphorus pollution risk assessment method based on migration and transformation mechanism according to claim 2, characterized in that: In step S20, the extracted physical, chemical and biological parameters are paired at the same location and time stamp to form a cross-media parameter group. The cross-media migration relationship of phosphorus is defined by the parameter group. The phosphorus speciation data in the chemical parameters is used as the core driving variable of the transformation relationship. The phosphorus speciation rate equation is constructed by combining the abundance of phosphorus cycling microorganisms in the biological parameters. Subsequently, the cross-medium parameter set is simultaneously input into the preset mechanism model and the preset output coefficient model, so that the spatial concentration field output by the preset mechanism model and the regional total load output by the preset output coefficient model are compared grid by grid on the same spatiotemporal grid to generate a dual-model residual matrix. The adjoint sensitivity algorithm is used to perform error inversion on the dual-model residual matrix, and the migration rate coefficient, conversion efficiency coefficient and correction coefficient of the preset mechanism model and the preset output coefficient model are updated synchronously at one time to realize the mutual constraint between micro-process and macro-accounting.

4. The phosphorus pollution risk assessment method based on migration and transformation mechanism according to claim 3, characterized in that: In step S30, the risk factor type screening is associated with the factor types in S20. Soluble exchangeable phosphorus, aluminum-bound phosphorus, and iron-bound phosphorus are the core chemical risk factors, phosphorus cycling microbial abundance is the core biological risk factor, and soil-water migration flux is the core physical risk factor. The output load standard is set according to the preset watershed ecological function zoning, which includes drinking water source protection areas and agricultural irrigation areas. Each level of standard meets the deviation rate threshold of the preset mechanism model and the preset output coefficient model and is dynamically adjusted according to soil type and water flow velocity. When processing the prediction results, the phosphorus pollution load accounting characteristics are weighted and coupled with the corresponding zoning load standard deviation and the proportion of high-value areas in the risk factor numerical distribution characteristics to form a comprehensive risk index.

5. The phosphorus pollution risk assessment method based on migration and transformation mechanism according to claim 4, characterized in that: In step S30, temporal characteristics include seasonal cycles, interannual variations, and extreme precipitation event windows, while spatial characteristics include upstream catchment areas, midstream main stream areas, downstream estuary areas, and key hydrological nodes such as reservoirs and dams. When predicting based on temporal characteristics, each migration model outputs phosphorus pollution load fluctuation curves at different time scales. The periods exceeding the standard are determined by combining the output load standards, and the corresponding risk factor numerical distributions are superimposed. When predicting based on spatial characteristics, the load thresholds for each sub-unit are first defined. After each migration model outputs the time series that meet the thresholds, the spatial correlation of risk factors in the sub-units at the time of compliance is extracted as a supplementary dimension to the prediction results.

6. The phosphorus pollution risk assessment method based on migration and transformation mechanism according to claim 5, characterized in that: In step S30, the phosphorus pollution load accounting characteristics are transformed into a load index based on the output load standard, and the risk factor numerical distribution characteristics are transformed into a three-level spatial risk level map of high, medium and low. Two types of features are used to generate a dynamic risk map through a spatiotemporal overlay algorithm. The map is marked with risk diffusion paths that meet the temporal characteristics and compliance time nodes that meet the spatial characteristics. Set an early warning threshold. When the product of the load index and the area ratio of high-risk zones exceeds the threshold, a risk warning will be automatically triggered. The warning information includes the most sensitive parameter type among the key coefficients.

7. The phosphorus pollution risk assessment method based on migration and transformation mechanism according to claim 6, characterized in that: In step S10, remote sensing data of the remote sensing area is collected by remote sensing equipment, vegetation data is extracted from the remote sensing images and processed into vegetation indices, and the biomass and phosphorus enrichment of the pre-set dominant plants are inverted by normalized vegetation indices; the phosphorus level of plant leaves is identified by hyperspectral remote sensing, and the density of benthic snails, waterbird activity and phosphorus input in feces are monitored simultaneously by infrared remote sensing. The obtained plant phosphorus enrichment and animal phosphorus input are included in the spatiotemporal matching database according to spatial coordinates and timestamps, and the local migration rate coefficient and conversion efficiency coefficient obtained by calibration with the pre-set mechanism model form a complementary verification of biological and chemical processes. In step S20, the relationship between soil plant retention and release is defined by the amount of phosphorus enrichment in plants, and the relationship between biological sedimentation and transformation in water bodies is constructed by animal population density and phosphorus input in feces. The preset mechanism model adds a plant phosphorus retention module and sets a phosphorus saturation absorption threshold. The preset output coefficient model incorporates animal phosphorus input coefficients. Error inversion simultaneously optimizes the plant phosphorus enrichment coefficient and animal phosphorus deposition coefficient, improving the accuracy of bio-mediated phosphorus cycle characterization and making up for the fragmented description of biological processes in traditional models.

8. The phosphorus pollution risk assessment method based on migration and transformation mechanism according to claim 7, characterized in that: In step S10, special collection is added at the preset sampling points. For human activities, the intensity of farmland fertilization, the frequency of phosphate mining, and the amount of domestic sewage discharge are collected and the corresponding interference intensity is set. For microorganisms, key soil layer samples for phosphorus loss are collected. The abundance of phosphate-solubilizing bacteria and polyphosphate-accumulating bacteria is determined by high-throughput sequencing and the phosphatase activity is detected by an enzyme activity analyzer. The interference intensity corresponding to human activities, the abundance of microorganisms, and the enzyme activity data are entered into the spatiotemporal matching database. In step S20, human activity data is correlated with a preset output coefficient model to correct the phosphorus loss coefficient under different disturbances. Microbial data is input into a preset mechanism model and an organophosphorus mineralization rate equation is constructed using microbial abundance and enzyme activity. A new human activity phosphorus release coefficient and microbial mineralization efficiency coefficient are added. Through error inversion and synchronous calibration with plant phosphorus enrichment coefficient and animal phosphorus deposition coefficient, a dual-dimensional correction of source strength identification and transformation simulation is achieved, enhancing the ability to characterize phosphorus migration and transformation in phosphate mining areas and farmland complex areas.

9. The phosphorus pollution risk assessment method based on migration and transformation mechanism according to claim 8, characterized in that: In step S10, remote sensing plant phosphorus enrichment data, animal phosphorus input data, and sampling point human activity disturbance intensity data, microbial abundance and enzyme activity data are paired according to the same location and time stamp to form a multi-dimensional parameter group of biological and human-made microorganisms, which is then incorporated into the spatiotemporal matching database and cross-validated with soil micro-characterization data. In step S20, the multi-dimensional parameter set is simultaneously input into the preset mechanism model and the preset output coefficient model. The preset mechanism model optimizes the adsorption-desorption-mineralization coupling parameters based on plant microbial data. The preset output coefficient model combines animal and anthropogenic data to correct the regional phosphorus load calculation formula. When generating the dual-model residual matrix, the residual of plant-microbial synergy and the residual of the interaction between anthropogenic mining and animal phosphorus deposition are included. The adjoint sensitivity algorithm is used to update the plant-microbial coupled phosphorus conversion coefficient and the anthropogenic-animal interaction phosphorus migration coefficient at one time.

10. A phosphorus pollution risk assessment system based on migration and transformation mechanisms, characterized in that, The method for assessing phosphorus pollution risk based on migration and transformation mechanism as described in any one of claims 1-9 was used.

Images