Data-free intelligent estimation method for ecological flow of mountainous area drainage basin based on multi-scale fusion

Abstract

The invention provides a data-free intelligent estimation method for ecological flow of a mountainous area drainage basin based on multi-scale fusion. The method comprises the following steps: acquiring basic data of a research area and preprocessing the basic data to obtain an effective data set; extracting a hydrological feature set and a process feature set according to the effective data set; generating a dynamic response unit based on the hydrological feature set and the process feature set, and performing multi-dimensional hydrological process decomposition to obtain a multi-scale feature matrix; performing ecological hydrological sequence matching on the multi-scale feature matrix to obtain a flow pattern set, and performing hierarchical uncertainty transfer analysis to generate an uncertainty index set; and constructing an integrated model according to the flow mode set and the uncertainty index set, and performing multi-constraint optimization to obtain an ecological flow estimation result. According to the method, the terrain complexity index and the hydrological similarity index are introduced, so that the problem of parameter uncertainty caused by spatial-temporal heterogeneity of rainfall in a mountainous area is effectively solved, and the accuracy of ecological flow estimation is improved.

Description

Technical Field

[0001] The invention relates to hydrological and hydrodynamic technology, in particular to an intelligent estimation method for ecological flow in a data-free mountainous area basin based on multi-scale fusion. Background Art

[0002] Intelligent estimation of ecological flow in data-free mountainous watersheds based on multi-scale fusion is of great significance for maintaining the health of river ecosystems and ensuring regional sustainable development. Mountainous watersheds are characterized by complex terrain, changeable climate, and strong heterogeneity of underlying surfaces. In addition, there are few observation sites, which poses a huge challenge to the accurate estimation of ecological flow. Accurately estimating ecological flow in mountainous watersheds can not only provide a scientific basis for ecological water allocation and river health assessment, but also provide decision-making support for integrated watershed management and ecological restoration, which is of great practical significance for balancing ecological protection and economic development.

[0003] At present, ecological flow estimation in mountainous watersheds mainly adopts hydrological methods, hydraulic methods and habitat methods. Hydrological methods determine ecological flow by setting fixed thresholds based on historical flow sequences, such as the Tennant method and the flow persistence curve method; hydraulic methods establish estimation models based on the relationship between hydraulic parameters and ecological needs, such as the wetted perimeter method and the R2CROSS method; habitat methods evaluate ecological flow by constructing habitat suitability models, such as the PHABSIM model. At the same time, researchers have also begun to try to introduce multi-source data such as remote sensing data and geographic information systems into ecological flow estimation, and use machine learning and other methods to improve estimation accuracy.

[0004] However, the existing research still has the following technical problems: the temporal and spatial heterogeneity of rainfall in mountainous areas is significant, and the existing methods are difficult to accurately characterize this heterogeneity in the parameter calibration process, resulting in large parameter uncertainty; in the process of multi-source data fusion, there are differences in the temporal and spatial resolutions of different data sources, and the traditional scale conversion method fails to effectively handle this scale difference, resulting in significant scale conversion errors; in terms of model integration, the existing methods mostly use static weight allocation strategies, which fail to fully adapt to the temporal and spatial heterogeneity of mountainous watersheds, affecting the reliability of the estimation results. The existence of these technical problems seriously restricts the accuracy and reliability of ecological flow estimation in mountainous watersheds without data. Summary of the invention

[0005] The purpose of the invention is to provide a method for intelligent estimation of ecological flow in data-free mountainous watersheds based on multi-scale fusion, in order to solve the above-mentioned problems existing in the conventional technology.

[0006] Technical solution: An intelligent estimation method for ecological flow in mountainous watersheds without data based on multi-scale fusion, including the following steps: Collect basic data of the study area and preprocess it to obtain effective data sets; basic data include DEM data, remote sensing image data, meteorological reanalysis data and land use data; Extract hydrological feature sets and process feature sets based on valid data sets; Generate dynamic response units based on hydrological feature sets and process feature sets, and perform multi-dimensional hydrological process decomposition to obtain a multi-scale feature matrix; The multi-scale characteristic matrix is ​​matched with eco-hydrological sequences to obtain the flow pattern set, and the hierarchical uncertainty transfer analysis is performed to generate the uncertainty indicator set; An integrated model was constructed based on the flow pattern set and uncertainty indicator set, and multi-constraint optimization was performed to obtain the ecological flow estimation results.

[0007] According to one aspect of the present application, the pre-processing comprises: Calculate the terrain complexity index value of the basic data to obtain the terrain decomposition factor; Based on the terrain decomposition factor, the basic data is locally decomposed and globally decomposed to generate a decomposition coefficient matrix; Calculate the hydrological similarity value according to the decomposition coefficient matrix to obtain the similarity index; Generate data quality indicators based on similarity indicators to obtain a valid data set.

[0008] According to one aspect of the present application, the step of extracting a hydrological feature set and a process feature set according to a valid data set includes: Extract the slope value, aspect value and elevation value of the grid unit, calculate the geomorphic position index and obtain the position index matrix; The geomorphic unit types are divided based on the position index matrix, and the characteristic parameters of each type of unit are extracted to generate the geomorphic feature vector; Constructing eigendecomposition basis functions based on topographic feature vectors to obtain a basis function set; Calculate the time response relationship between rainfall data and runoff data to obtain the time lag coefficient matrix; According to the basis function set and the time-delay coefficient matrix, the effective data set is subjected to feature extraction to obtain the feature component set; The characteristic component sets are combined at multiple scales to generate hydrological characteristic sets and process characteristic sets.

[0009] According to one aspect of the present application, the step of generating a dynamic response unit includes: Calculate the terrain gradient value of each point in the watershed to obtain the gradient matrix; Construct spatial neighborhood relationships based on the gradient matrix and generate a neighborhood matrix; Calculate the confluence time of each area and obtain the time matrix; A dynamic partitioning standard is constructed based on the neighborhood matrix and the time matrix to obtain the partitioning criterion; Perform adaptive partitioning according to the partitioning criteria to obtain the initial partitioning result; The initial partitioning results are dynamically adjusted according to seasonal changes to obtain dynamic response units.

[0010] According to one aspect of the present application, the steps of performing multidimensional hydrological process decomposition include: Extracting the terrain features of each unit to obtain a terrain feature set; Calculate the hydrological response parameters of each unit to obtain a response parameter set; Construct the characteristic expression of unit scale to obtain the characteristic function set; Perform multi-dimensional decomposition based on the characteristic function set to obtain an initial decomposition set; Optimize the decomposition coefficients of the initial decomposition set to generate a decomposition coefficient set; The energy value of each decomposition scale is calculated based on the decomposition coefficient set to obtain a multi-scale feature matrix.

[0011] According to one aspect of the present application, the step of obtaining a flow pattern set by matching an eco-hydrological sequence includes: Extract the water volume change sequence to obtain the water volume sequence; Calculate the characteristics of ecological sensitive periods to obtain the sensitive period matrix; Analyze ecological water demand to obtain the demand matrix; The response intensity is calculated based on the water quantity series and the sensitive period matrix to obtain the intensity matrix; Combining the demand matrix and the intensity matrix to construct response characteristics and generate a response feature set; Pattern recognition is performed based on the response feature set to generate a traffic pattern set.

[0012] According to one aspect of the present application, the steps of performing hierarchical uncertainty transfer analysis to generate an uncertainty indicator set include: Construct the hydrological process transmission chain network to obtain the network structure; Calculate the parameter sensitivity index of each node to generate a sensitivity matrix; Analyze the uncertainty transfer path according to the network structure and sensitivity matrix to obtain the transfer path set; Evaluate the coupling effect of multi-source uncertainties to obtain a set of coupling coefficients; The comprehensive uncertainty is calculated based on the transfer path set and the coupling coefficient set to generate an uncertainty index set.

[0013] According to one aspect of the present application, the step of constructing an integrated model includes: Extract terrain feature parameters to obtain terrain feature matrix; Calculate the hydrological response index to obtain the response index matrix; Analyze land use characteristics to obtain land characteristic matrix; Assess climate variability to obtain a climate indicator matrix; Model selection criteria were constructed based on terrain characteristics matrix, response indicator matrix, land characteristics matrix and climate indicator matrix; Calculate the spatiotemporal heterogeneity index of each region to obtain the heterogeneity matrix; Determine the model weights according to the model selection criteria and the heterogeneity matrix to generate a weight vector; Evaluate the reliability index of each model to obtain the reliability matrix; The integrated model is constructed using the weight vector and reliability matrix.

[0014] According to one aspect of the present application, a multi-constraint optimization is performed to obtain an ecological flow estimation result. The steps include: Construct the seasonal variation function to obtain the variation function; Design the ecological demand constraint system to obtain the constraint condition set; Construct an objective function that takes seasonal changes into account; Perform optimization and solution based on the constraint condition set and the objective function to obtain the initial result; Analyze the response characteristics of the ecosystem to the initial results to obtain response indicators; The initial results are revised based on the response indicators to generate ecological flow estimates.

[0015] According to one aspect of the present application, ecohydrological sequence matching further includes: Extract river network structure data to obtain river network matrix; Calculate the physical connectivity between river sections to obtain the physical connectivity matrix; Analyze the hydrological connectivity to obtain the hydrological connectivity matrix; Evaluate the ecological corridor function to obtain the corridor index matrix; The physical connectivity matrix, hydrological connectivity matrix and corridor index matrix were integrated to generate the connectivity matrix; The connectivity matrix, together with the response feature set, is used for pattern recognition to generate a set of traffic patterns.

[0016] According to one aspect of the present application, the step of analyzing the uncertainty transfer path includes: Design a perturbation experiment plan to obtain a perturbation plan; Perform parameter perturbation test to obtain test results; Calculate the response sensitivity to obtain the sensitivity matrix; Analyze the importance of parameters to obtain the importance vector; Identify the sources of uncertainty and obtain the source matrix; Analyze the coupling action mode to obtain the action mode matrix; Generate a set of coupling coefficients based on the action mode matrix.

[0017] Beneficial effect: By introducing the terrain complexity index and hydrological similarity index, the parameter uncertainty problem caused by the spatiotemporal heterogeneity of rainfall in mountainous areas was effectively solved, and the accuracy of ecological flow estimation was improved. BRIEF DESCRIPTION OF THE DRAWINGS

[0018] Figure 1 It is a flow chart of the present invention.

[0019] Figure 2 It is a flow chart of the preprocessing of the present invention.

[0020] Figure 3 It is a flow chart of the present invention for extracting a hydrological feature set and a process feature set based on a valid data set.

[0021] Figure 4 It is a flow chart of generating a dynamic response unit according to the present invention.

[0022] Figure 5 It is a flow chart of the multi-dimensional hydrological process decomposition of the present invention. DETAILED DESCRIPTION

[0023] like Figure 1 As shown in the figure, the intelligent estimation method of ecological flow in mountainous watershed without data based on multi-scale fusion includes: Step S1: Receive DEM data, remote sensing image data, meteorological reanalysis data, and land use data, use a multi-scale adaptive decomposition method to process various input data, generate data quality indicators and valid data sets; based on the valid data sets, perform adaptive extraction of hydrological process characteristics to obtain a hydrological feature set and a process feature set.

[0024] Step S2: receiving a hydrological feature set and a process feature set, processing them using a terrain-guided adaptive partitioning method, and generating a dynamic response unit; performing a multi-dimensional hydrological process decomposition based on the dynamic response unit to obtain a multi-scale feature matrix.

[0025] Step S3: Receive the multi-scale feature matrix, process it using the eco-hydrological sequence matching algorithm, and obtain the flow pattern set; perform hierarchical uncertainty transfer analysis based on the flow pattern set to generate an uncertainty indicator set.

[0026] Step S4: Receive the flow pattern set and the uncertainty indicator set, and use the adaptive dynamic integration method to build the integrated model E; perform multi-constraint optimization based on the integrated model to obtain the ecological flow estimation results and reliability evaluation indicators.

[0027] In addition to the problems in the background technology, the study found that the ecosystem exhibits obvious nonlinear response characteristics to flow changes. Existing methods mostly use linear or simple nonlinear models to describe them, which makes it difficult to accurately characterize this complex response relationship; in addition, under extreme climatic conditions, ecological flow shows unique dynamic changes, but existing research mainly focuses on flow characteristics under conventional climatic conditions, and lacks in-depth analysis of the dynamic characteristics of ecological flow under extreme conditions; finally, current uncertainty analysis methods often separate each link and lack a systematic analysis of the uncertainty transmission mechanism in the hydrological process chain, making it difficult to quantify the coupling effect of multi-source uncertainty.

[0028] Through the technical route of multi-scale fusion, the intelligent estimation of ecological flow in mountainous watersheds without data is realized. The scheme effectively solves the technical problems of scarce data, complex terrain, and variable process response in mountainous watersheds through innovative methods such as multi-scale adaptive decomposition, terrain-guided dynamic partitioning, eco-hydrological sequence matching, and adaptive dynamic integration. In particular, when processing multi-source heterogeneous data, an adaptive decomposition strategy based on terrain features is adopted to improve data quality and reliability; in the process of model construction, a refined description of the hydrological process in mountainous watersheds is achieved through dynamic response unit division and multi-dimensional feature extraction; in the link of ecological flow estimation, the combination of ecological demand constraints and multi-objective optimization ensures the scientificity and reliability of the estimation results. The overall scheme has strong universality and practicality, and can provide strong technical support for ecological flow management in mountainous watersheds.

[0029] According to one aspect of the present application, step S1 specifically comprises: Step S11: receiving DEM data, remote sensing image data, meteorological reanalysis data, and land use data, calculating the terrain complexity index value of each data, and obtaining the terrain decomposition factor; performing local and global decomposition on various input data based on the terrain decomposition factor, and generating a decomposition coefficient matrix; calculating the hydrological similarity value between each data point and adjacent data points according to the decomposition coefficient matrix, and obtaining a similarity index; using the similarity index to identify outliers, and generating a data quality index; screening the original input data based on the data quality index, and obtaining a valid data set.

[0030] Step S12: Receive a valid data set, extract characteristic parameters of each watershed geomorphic unit, and construct a geomorphic characteristic vector; generate characteristic decomposition basis functions according to the geomorphic characteristic vector to obtain a basis function set; calculate the time response relationship between rainfall data and runoff data to obtain a time lag coefficient matrix; use the basis function set and the time lag coefficient matrix to extract features of the valid data set to obtain a feature component set; perform multi-scale combination of the feature component set to generate a hydrological feature set and a process feature set.

[0031] The input data is processed by a multi-scale adaptive decomposition method, which realizes the efficient fusion and quality control of multi-source heterogeneous data under complex terrain conditions in mountainous areas. Specifically, the spatial variability of terrain characteristics is quantified by calculating the terrain complexity index. The terrain decomposition factor constructed based on this can adaptively adjust the scale and intensity of data decomposition, so that a more detailed decomposition method is used in complex terrain areas such as steep mountains and canyons, while a relatively rough decomposition is used in areas with flat terrain, thereby ensuring the matching degree between the decomposition results and the actual terrain characteristics. At the same time, by calculating the hydrological similarity between data points, an outlier identification standard based on physical mechanisms is established, which can effectively screen out abnormal data caused by factors such as terrain and meteorology, and improve the data quality of subsequent analysis. This data preprocessing method based on multi-scale decomposition not only improves the reliability of the data, but also the hydrological feature set and process feature set obtained by feature extraction can comprehensively characterize the hydrological process characteristics of mountainous watersheds, providing a high-quality data basis for subsequent ecological flow estimation.

[0032] According to one aspect of the present application, step S2 is specifically: Step S21: Receive the hydrological feature set and the process feature set, calculate the terrain gradient value of each point in the watershed, and obtain the gradient matrix; construct the spatial neighborhood relationship according to the gradient matrix to generate the neighborhood matrix; calculate the confluence time of each area to obtain the time matrix; construct the dynamic partitioning standard based on the neighborhood matrix and the time matrix to obtain the partitioning criterion; use the partitioning criterion to adaptively partition the study area to generate the initial partitioning result; dynamically adjust the initial partitioning result according to seasonal changes to obtain the dynamic response unit.

[0033] Step S22: Receive dynamic response units, extract watershed characteristic parameters of each unit, and construct a characteristic function set; perform multi-dimensional decomposition on each response unit according to the characteristic function set to obtain a decomposition coefficient set; calculate the energy distribution of each component in the decomposition coefficient set to obtain an energy matrix; identify the dominant mode in the energy matrix and generate a mode matrix; calculate the information flow between each scale to obtain a flow index; perform feature reorganization according to the mode matrix and the flow index to generate a multi-scale feature matrix.

[0034] The terrain-guided adaptive partitioning method is used to achieve a refined description and analysis of complex hydrological processes in mountainous watersheds. This method first constructs spatial neighborhood relationships based on terrain gradients, and establishes partitioning standards that take into account hydrological connectivity by calculating confluence time, so that the partitioning results can accurately reflect the hydrological response characteristics of mountainous watersheds. Especially in mountainous areas with significant seasonal changes, by dynamically adjusting the partition boundaries, the changing characteristics of the watershed hydrological process in different seasons can be effectively captured. This dynamic response unit division method breaks through the limitations of traditional static partitioning methods and can better adapt to the spatiotemporal heterogeneity of hydrological processes in mountainous watersheds. At the same time, through multidimensional hydrological process decomposition, complex hydrological processes are decomposed into characteristic components of different scales, and the dominant mode is identified through energy distribution analysis, realizing the multi-scale feature extraction of hydrological processes in mountainous watersheds, providing a scientific basis for accurately estimating ecological flow.

[0035] According to one aspect of the present application, step S3 is specifically: Step S31: receiving a multi-scale feature matrix, extracting the ecological response characteristics of each time series, and obtaining a response feature set; calculating the connectivity index of the water system in the basin, and generating a connectivity matrix; constructing a similarity metric based on the response feature set and the connectivity matrix, and obtaining a similarity index; performing multi-time scale decomposition on the data sequence, and obtaining a time series set; performing pattern matching on the time series set using the similarity index, and generating a flow pattern set.

[0036] Step S32: Receive the flow pattern set, construct the hydrological process transmission chain network, and obtain the network structure; calculate the parameter sensitivity index of each node and generate a sensitivity matrix; analyze the uncertainty transmission path according to the network structure and the sensitivity matrix to obtain the transmission path set; evaluate the coupling effect of multi-source uncertainty to obtain the coupling coefficient set; calculate the comprehensive uncertainty based on the transmission path set and the coupling coefficient set to generate the uncertainty index set.

[0037] Through the eco-hydrological sequence matching algorithm and hierarchical uncertainty transmission analysis, the accurate identification and uncertainty quantification of ecological flow characteristics in mountainous watersheds are achieved. This method first constructs a similarity metric based on ecological response characteristics and water system connectivity. Through multi-time scale decomposition and pattern matching, it can accurately identify the flow pattern characteristics under different hydrological scenarios. Especially in mountainous watersheds where data is scarce, this sequence matching-based method can make full use of limited observational data to extract representative flow patterns. At the same time, by constructing a hydrological process transmission chain network, analyzing the transmission rules of uncertainty between different levels, and evaluating the coupling effect of multi-source uncertainty, the uncertainty of ecological flow estimation results is systematically quantified, providing a reliable uncertainty assessment basis for decision-making.

[0038] According to one aspect of the present application, step S4 is specifically: Step S41: receiving a flow pattern set and an uncertainty indicator set, extracting basin characteristic parameters, and generating a characteristic indicator set; establishing a model selection standard based on the characteristic indicator set to obtain a selection criterion; calculating the spatiotemporal heterogeneity index of each region to obtain a heterogeneity matrix; determining the model weight based on the selection criterion and the heterogeneity matrix to generate a weight vector; evaluating the reliability index of each model to obtain a reliability matrix; constructing an integrated model using the weight vector and the reliability matrix to generate an integrated model.

[0039] Step S42: Receive the integrated model, establish an ecological demand constraint system, and obtain a set of constraint conditions; construct an objective function that takes seasonal changes into account and generate an objective function; perform optimization and solution based on the constraint set and the objective function to obtain an initial result; analyze the response characteristics of the ecosystem to the initial result and obtain a response index; modify the initial result based on the response index and generate an ecological flow estimation result; evaluate the reliability of the ecological flow estimation result and obtain a reliability evaluation index.

[0040] Through adaptive dynamic integration method and multi-constraint optimization, accurate estimation and reliability assessment of ecological flow in mountainous watersheds are achieved. This method first constructs model selection criteria based on watershed characteristics, determines model weights considering temporal and spatial heterogeneity, and establishes an integrated model that adapts to the characteristics of mountainous watersheds. Especially when considering seasonal changes, by establishing an ecological demand constraint system and multi-objective optimization, it is possible to ensure the scientific rationality of the estimation results while meeting the needs of the ecosystem. At the same time, by analyzing the response characteristics of the ecosystem, the estimation results are corrected, and a systematic reliability assessment system is established, which realizes dynamic optimization and quality control of the estimation results, providing reliable technical support for ecological flow management in mountainous watersheds.

[0041] According to one aspect of the present application, step S12 is specifically: Step S121: Receive a valid data set, extract the slope value of the grid cell, and obtain a slope matrix; calculate the aspect value of the grid cell to obtain an aspect matrix; extract the elevation value of the grid cell to obtain an elevation matrix; calculate the landform position index based on the slope matrix, the aspect matrix and the elevation matrix to obtain a position index matrix; divide the landform unit types based on the position index matrix to obtain a unit type matrix; extract parameters such as the area ratio, average elevation and slope of each type of unit based on the unit type matrix to generate a landform feature vector.

[0042] Step S122: Receive the geomorphic feature vector, construct the spatial distribution function of the geomorphic unit, and obtain the spatial function set; calculate the confluence path length of each geomorphic unit to obtain the path length matrix; construct the basis function expression according to the spatial function set and the path length matrix to generate the basis function set; extract the rainfall data sequence from the valid data set to obtain the rainfall sequence; extract the runoff data sequence from the valid data set to obtain the runoff sequence; calculate the mutual correlation coefficient between the rainfall sequence and the runoff sequence to generate the correlation coefficient matrix; identify the peak response time based on the correlation coefficient matrix to obtain the time lag coefficient matrix R.

[0043] Step S123: receiving a basis function set, a time-lag coefficient matrix and a valid data set, constructing a projection matrix for feature extraction, and obtaining a projection matrix; transforming the valid data set using the projection matrix to obtain a transformation coefficient set; adjusting the transformation coefficient set in time scale according to the time-lag coefficient matrix to obtain an adjustment coefficient set; calculating the principal components of the adjustment coefficient set to generate a feature component set.

[0044] Step S124: receiving a feature component set, constructing a multi-scale decomposition window sequence, and obtaining a window sequence; performing sliding decomposition on the feature component set according to the window sequence to obtain a decomposition coefficient set; calculating the energy distribution of the decomposition coefficient set at different scales to obtain an energy distribution matrix; selecting a scale threshold for feature recombination according to the energy distribution matrix to obtain a threshold vector; recombining the decomposition coefficient set based on the threshold vector to generate a hydrological feature set; extracting process feature parameters according to the hydrological feature set to obtain a process feature set.

[0045] Through the feature extraction and multi-scale combination method based on geomorphic units, the systematic extraction and characterization of the hydrological characteristics of mountainous watersheds are realized. This method first calculates the geomorphic position index based on terrain factors such as slope, aspect and elevation, constructs the characteristic basis function through spatial distribution function and confluence path analysis, and establishes a feature extraction framework that adapts to the terrain characteristics of mountainous areas. In particular, when dealing with the rainfall-runoff response relationship, the peak response time is identified by calculating the mutual correlation coefficient, and the hydrological response characteristics of mountainous watersheds are accurately characterized. At the same time, multi-scale sliding decomposition and energy distribution analysis are used to realize the multi-scale reorganization of hydrological characteristics, which not only retains the important feature information at different scales, but also improves the accuracy and representativeness of feature expression. This feature extraction method fully considers the terrain characteristics and hydrological process characteristics of mountainous watersheds, and provides reliable feature support for subsequent ecological flow estimation.

[0046] According to one aspect of the present application, step S21 is specifically: Step S211: Receive the hydrological feature set and the process feature set, extract the terrain elevation data, and obtain the elevation matrix; calculate the height difference between adjacent grids in the elevation matrix to obtain the elevation difference matrix; calculate the slope value according to the elevation difference matrix and the grid spatial resolution to obtain the slope matrix; calculate the terrain gradient by combining the slope matrix and the flow direction information to generate the gradient matrix.

[0047] Step S212: Receive the gradient matrix, define the spatial search radius, and obtain the search parameters; identify the neighborhood units of each grid according to the search parameters to obtain a neighborhood set; calculate the spatial connection strength between the units in the neighborhood set to obtain a connection strength matrix; construct a neighborhood relationship network based on the connection strength matrix to generate a neighborhood matrix.

[0048] Step S213: Receive the neighborhood matrix, extract the river network structure information, and obtain the river network matrix; calculate the flow velocity parameters of each river section in the river network matrix to obtain the flow velocity matrix; calculate the confluence path by combining the gradient matrix and the flow velocity matrix to obtain the path matrix; calculate the confluence time based on the path matrix to generate the time matrix.

[0049] Step S214: receiving the neighborhood matrix and the time matrix, constructing a spatial constraint function, and obtaining a constraint function; calculating the similarity index between units according to the constraint function, and obtaining a similarity matrix; setting a partition threshold standard, and obtaining a threshold vector; and constructing a partition criterion by combining the similarity matrix and the threshold vector to generate a partition criterion.

[0050] Step S215: receiving partition criteria, initializing partition identifiers, and obtaining an identifier matrix; C Iteratively optimize the identification matrix to obtain the optimized matrix; test the spatial continuity of the optimized matrix to obtain the continuity index; modify the partition result based on the continuity index to generate the initial partition result.

[0051] Step S216: Receive the initial partitioning results, extract seasonal hydrological characteristics, and obtain a seasonal characteristic matrix; construct a seasonal adjustment function to obtain an adjustment function; calculate the dynamic change range of the partition boundary according to the adjustment function to obtain a change range matrix; adjust the initial partitioning results based on the change range matrix to generate a dynamic response unit.

[0052] Through the terrain-guided adaptive partitioning method, dynamic partitioning and spatial unit division of mountainous watersheds are realized. This method first calculates the terrain gradient based on elevation and slope, establishes a neighborhood relationship network through spatial search and connection strength analysis, and constructs a spatial partitioning framework that takes terrain characteristics into consideration. In particular, when calculating the confluence time, the river network structure and flow velocity parameters are combined to accurately reflect the hydrological connectivity characteristics of mountainous watersheds. At the same time, by constructing spatial constraint functions and similarity indicators, the partitioning results are optimized and adjusted, and the partition boundaries are dynamically adjusted through seasonal characteristic analysis, so that the partitioning results can better adapt to the spatiotemporal variation characteristics of the hydrological process in mountainous watersheds. This dynamic partitioning method breaks through the limitations of traditional static partitioning and provides a more accurate basis for spatial unit division for the analysis of hydrological processes in mountainous watersheds.

[0053] According to one aspect of the present application, step S22 is specifically: Step S221: receiving dynamic response units, extracting terrain features of each unit, and obtaining a terrain feature set; calculating hydrological response parameters of each unit, and obtaining a response parameter set; constructing a characteristic expression at the unit scale, and obtaining an expression matrix; generating a characteristic function according to the expression matrix, and obtaining a characteristic function set.

[0054] Step S222: receiving a characteristic function set, constructing a multidimensional decomposition basis function, and obtaining a basis function matrix; designing a decomposition scale sequence, and obtaining a scale sequence; performing multidimensional decomposition according to the basis function matrix and the scale sequence, and obtaining an initial decomposition set; optimizing the decomposition coefficients of the initial decomposition set, and generating a decomposition coefficient set.

[0055] Step S223: receiving a decomposition coefficient set, calculating the energy value of each decomposition scale, and obtaining an energy value matrix; constructing an energy distribution function, and obtaining a distribution function; using the distribution function to analyze energy distribution characteristics, and generating an energy matrix.

[0056] Step S224: receiving the energy matrix, setting the pattern recognition criteria, and obtaining the recognition criteria; extracting the dominant pattern features according to the recognition criteria, and obtaining the feature matrix; analyzing the spatiotemporal distribution of the dominant pattern, and obtaining the distribution matrix; integrating the feature matrix and the distribution matrix, and generating the pattern matrix.

[0057] Step S225: receiving the pattern matrix, constructing the scale conversion function, and obtaining the conversion function; calculating the amount of information transmission between adjacent scales, and obtaining the transmission matrix; evaluating the intensity of information flow, and obtaining the intensity index; and combining the transmission matrix and the intensity index to generate the flow index.

[0058] Step S226: receiving the pattern matrix and the traffic index, determining the feature recombination weight, and obtaining the weight vector; performing feature recombination calculation according to the weight vector, and obtaining the recombination matrix; evaluating the representativeness of the recombination result, and obtaining the representative index; optimizing the recombination matrix based on the representative index, and generating a multi-scale feature matrix.

[0059] Through the multidimensional hydrological process decomposition and feature reorganization method, the multi-scale feature extraction of complex hydrological processes in mountainous watersheds is realized. This method first constructs characteristic functions based on terrain characteristics and hydrological response parameters, decomposes the process through multidimensional decomposition basis functions and decomposition scale sequences, and establishes a feature extraction framework that adapts to the characteristics of mountainous watersheds. In particular, in the process of energy distribution analysis, the key features of the hydrological process are accurately captured by identifying the dominant mode and analyzing the spatiotemporal distribution characteristics. At the same time, by calculating the information flow between scales and constructing feature reorganization weights, the optimized reorganization of the decomposed features is achieved, which not only ensures the representativeness of the reorganization results, but also improves the accuracy of feature expression. This multidimensional decomposition method can effectively handle the multi-scale characteristics of hydrological processes in mountainous watersheds and provide reliable feature support for ecological flow estimation.

[0060] According to one aspect of the present application, step S31 is specifically: Step S311: Receive a multi-scale feature matrix, extract a water volume change sequence, and obtain a water volume sequence; calculate ecological sensitive period characteristics to obtain a sensitive period matrix; analyze ecological water demand to obtain a demand matrix; calculate the response intensity based on the water volume sequence and the sensitive period matrix to obtain an intensity matrix; combine the demand matrix and the intensity matrix to construct response characteristics and generate a response feature set.

[0061] Step S312: Receive the multi-scale feature matrix, extract the river network structure data, and obtain the river network matrix; calculate the physical connectivity between river sections to obtain the physical connectivity matrix; analyze the hydrological connectivity status to obtain the hydrological connectivity matrix; evaluate the ecological corridor function to obtain the corridor index matrix; integrate the physical connectivity matrix, the hydrological connectivity matrix and the corridor index matrix to generate the connectivity matrix C .

[0062] Step S313: receiving the response feature set and the connectivity matrix, constructing a feature weight function, and obtaining a weight function; calculating the correlation between features, and obtaining a correlation matrix; constructing a distance function based on the weight function and the correlation matrix, and obtaining a distance function; designing a similarity calculation criterion, and obtaining a calculation criterion; generating a similarity index S according to the distance function and the calculation criterion.

[0063] Step S314: receiving a multi-scale feature matrix, determining a time window sequence, and obtaining a window sequence; constructing a time scale decomposition function, and obtaining a decomposition function; performing sequence decomposition according to the window sequence and the decomposition function, and obtaining a decomposition result; and reorganizing and optimizing the decomposition result to generate a time series set.

[0064] Step S315: Receive similarity index and time series set, construct pattern matching criteria, and obtain matching criteria; calculate dynamic similarity between sequences to obtain dynamic similarity matrix; identify key matching points to obtain matching point set; perform pattern recognition based on matching criteria and matching point set to generate traffic pattern set.

[0065] Through the eco-hydrological sequence matching and response characteristic analysis method, the accurate identification of ecological flow patterns in mountainous watersheds is achieved. This method first analyzes the ecological response characteristics based on water volume changes and ecological sensitive period characteristics, establishes hydrological connectivity through river network connectivity and ecological corridor function assessment, and constructs a pattern recognition framework that takes into account ecological needs. In particular, in the similarity calculation process, the dynamic similarity relationship between flow sequences is accurately characterized by combining feature weights and correlation analysis. At the same time, through time scale decomposition and pattern matching, accurate identification of flow patterns is achieved, which not only ensures the reliability of identification results, but also improves the accuracy of pattern expression. This pattern recognition method based on ecological response fully considers the ecological characteristics of mountainous watersheds and provides a scientific model support for ecological flow estimation.

[0066] According to one aspect of the present application, step S32 is specifically: Step S321: Receive a traffic pattern set, extract key process nodes, and obtain a node set; analyze the connection relationship between nodes to obtain a connection matrix; construct a process transfer function to obtain a transfer function; establish a network topology based on the node set and the connection matrix to obtain a topology matrix; construct a transfer chain network based on the transfer function and the topology matrix to generate a network structure.

[0067] Step S322: Receive the traffic pattern set and network structure, design a disturbance experiment plan, and obtain a disturbance plan; perform a parameter disturbance test to obtain test results; calculate the response sensitivity to obtain a sensitivity matrix; analyze the parameter importance to obtain an importance vector; integrate the sensitivity matrix and the importance vector to generate a sensitivity matrix.

[0068] Step S323: Receive the network structure and sensitivity matrix, construct the transfer path identification criteria, and obtain the identification criteria; analyze the uncertainty propagation characteristics to obtain the propagation characteristic matrix; identify the key transfer path to obtain the key path set; evaluate the path importance to obtain the importance matrix; integrate the key path set and the importance matrix to generate the transfer path set.

[0069] Step S324: Receive the transfer path set, identify the uncertainty source, and obtain the source matrix; construct the coupling effect evaluation function to obtain the evaluation function; calculate the uncertainty interaction intensity to obtain the interaction matrix; analyze the coupling action mode to obtain the action mode matrix; generate the coupling coefficient set based on the interaction matrix and the action mode matrix C .

[0070] Step S325: Receive the transmission path set and the coupling coefficient set, construct a comprehensive evaluation function, and obtain an evaluation function; calculate the uncertainty contribution of each path to obtain a contribution matrix; evaluate the overall uncertainty level to obtain a level index; analyze the spatiotemporal distribution of uncertainty to obtain a distribution matrix; integrate the level index and the distribution matrix to generate an uncertainty index set.

[0071] Through uncertainty transmission analysis and coupling effect evaluation methods, the uncertainty of ecological flow estimation in mountainous watersheds is systematically quantified. This method first constructs a transmission chain network based on key nodes of the process, identifies key influencing factors through parameter perturbation tests and sensitivity analysis, and establishes an uncertainty analysis framework that adapts to the characteristics of mountainous watersheds. In particular, when analyzing the characteristics of uncertainty propagation, the transmission law of uncertainty is accurately characterized by identifying key transmission paths and evaluating the importance of paths. At the same time, by analyzing the coupling effects of multi-source uncertainties and constructing a comprehensive evaluation function, the systematic quantification of uncertainty is achieved, which not only ensures the comprehensiveness of the evaluation results, but also improves the accuracy of uncertainty characterization. This uncertainty analysis method can effectively evaluate the reliability of ecological flow estimation in mountainous watersheds and provide important support for decision-making.

[0072] According to one aspect of the present application, step S41 is specifically: Step S411: Receive a flow pattern set and an uncertainty index set, extract terrain characteristic parameters, and obtain a terrain characteristic matrix; calculate hydrological response indicators to obtain a response index matrix; analyze land use characteristics to obtain a land characteristic matrix; evaluate climate variability to obtain a climate index matrix; integrate the terrain characteristic matrix, the response index matrix, the land characteristic matrix, and the climate index matrix to generate a characteristic index set.

[0073] Step S412: Receive a feature indicator set, construct a model adaptability evaluation function, and obtain an evaluation function; calculate feature indicator weights to obtain an indicator weight matrix; design a model selection threshold to obtain a threshold vector; construct a selection rule based on the evaluation function and the indicator weight matrix to obtain a selection rule; generate a selection criterion based on the selection rule and the threshold vector.

[0074] Step S413: receiving a characteristic indicator set, dividing the spatiotemporal calculation units, and obtaining a calculation unit set; analyzing the spatial heterogeneity characteristics to obtain a spatial heterogeneity matrix; calculating the temporal variability index to obtain a temporal variability matrix; evaluating the scale dependence to obtain a scale indicator matrix; integrating the spatial heterogeneity matrix, the temporal variability matrix, and the scale indicator matrix to generate a heterogeneity matrix.

[0075] Step S414: Receive the selection criteria and heterogeneity matrix, construct a weight distribution function, and obtain a distribution function; calculate the model adaptability score to obtain a score matrix; optimize the weight distribution scheme to obtain an optimization scheme; calculate the weight according to the distribution function DF and the optimization scheme to generate a weight vector. Step S415: Receive a traffic pattern set, design a verification test plan, and obtain a test plan; perform model verification calculations to obtain verification results; evaluate prediction accuracy to obtain an accuracy matrix; analyze stability indicators to obtain stability indicators; integrate the accuracy matrix and stability indicators to generate a reliability matrix.

[0076] Step S416: receiving the weight vector and the reliability matrix, constructing a model integration framework, and obtaining an integrated framework; designing a model combination strategy, and obtaining a combination strategy; optimizing the integration parameters, and obtaining a parameter matrix; constructing an integrated model based on the integrated framework and the combination strategy, and generating an integrated model.

[0077] Through adaptive dynamic integration and model selection methods, the optimal integration of ecological flow estimation models in mountainous watersheds is achieved. This method first constructs a characteristic indicator set based on factors such as terrain characteristics, hydrological response, land use and climate variability, establishes an integration framework through model adaptability evaluation and weight optimization, and builds a model integration system that adapts to the characteristics of mountainous watersheds. Especially when considering spatiotemporal heterogeneity, the applicability characteristics of the model are accurately characterized by analyzing spatial heterogeneity, temporal variability and scale dependence. At the same time, through model verification and reliability evaluation, the optimization adjustment of the integrated model is achieved, which not only ensures the prediction accuracy of the model, but also improves the reliability of the estimation results. This adaptive integration method can make full use of the advantages of multiple models and provide a reliable model support for ecological flow estimation in mountainous watersheds.

[0078] According to one aspect of the present application, step S42 is specifically: Step S421: Receive the integrated model, extract ecological demand characteristics, and obtain a demand characteristic matrix; construct water constraint conditions to obtain a water constraint matrix; design ecological protection goals to obtain a protection target matrix; analyze the system carrying capacity to obtain a carrying capacity index; integrate the water constraint matrix, the protection target matrix and the carrying capacity index to generate a constraint condition set.

[0079] Step S422: receiving a set of constraints, constructing a seasonal variation function, and obtaining a variation function; designing a target optimization criterion, and obtaining an optimization criterion; establishing a multi-objective trade-off function, and obtaining a trade-off function; constructing a target function based on the variation function and the trade-off function, and generating a target function.

[0080] Step S423: receiving a set of constraints and an objective function, designing optimization algorithm parameters, and obtaining a parameter set; performing optimization iterative calculations, and obtaining iterative results; analyzing convergence characteristics, and obtaining convergence indicators; selecting the optimal solution according to the convergence indicators, and generating an initial result.

[0081] Step S424: receiving the initial results, constructing an ecological response evaluation function, and obtaining an evaluation function; calculating the change of ecological indicators, and obtaining a change matrix; analyzing the response sensitivity, and obtaining a sensitivity index; evaluating the response characteristics, and obtaining a characteristic matrix; integrating the sensitivity index and the characteristic matrix, and generating a response index.

[0082] Step S425: Receive the initial results and response indicators, construct a correction function model, and obtain the correction function; calculate the adjustment parameters to obtain the adjustment matrix; perform result correction to obtain the correction result; evaluate the correction effect to obtain the effect indicator; determine the final result based on the effect indicator and generate the ecological flow estimation result.

[0083] Step S426: Receive the ecological flow estimation results, construct an evaluation index system, and obtain the index system; calculate the reliability index to obtain the reliability matrix; analyze the uncertainty level to obtain the uncertainty index; evaluate the adaptability characteristics to obtain the adaptability matrix; integrate the reliability matrix, uncertainty index and adaptability matrix to generate the reliability evaluation index.

[0084] Through multi-constraint optimization and ecological response evaluation methods, dynamic optimization of ecological flow estimation results in mountainous watersheds is achieved. This method first constructs constraints based on ecological demand characteristics, establishes an optimization framework through seasonal variation functions and multi-objective trade-offs, and builds an optimization system that takes into account ecosystem needs. In particular, during the optimization process, the optimal solution is accurately identified through iterative calculation and convergence analysis, and the results are corrected through ecological response evaluation. At the same time, by constructing an evaluation index system and analyzing reliability characteristics, a systematic evaluation of the estimation results is achieved, which not only ensures the scientific nature of the results, but also improves the comprehensiveness of the evaluation. This multi-constraint optimization method can effectively balance ecological needs and system constraints, and provides a scientific basis for ecological flow management in mountainous watersheds.

[0085] Case 1: Step S1: Multi-source heterogeneous data preprocessing and feature extraction The collected data include: DEM data, 30m resolution grid data, including elevation value H(i, j); remote sensing image data, 10m resolution Sentinel-2 image, including multispectral band value B(i, j, k); meteorological reanalysis data, including daily data such as rainfall P(t) and temperature T(t); land use data, 30m resolution grid data, including land use type code L(i, j); Step S11: Specific process of data quality assessment For each grid point (i, j), extract the elevation value within the 3×3 window; calculate the terrain complexity index T C (i, j):T C (i,j)=√[(∑(H(i,j)-H_mean) 2 ) / n]; where H_mean is the average elevation in the window, and n is the number of grids in the window; the terrain decomposition factor T is generated in the form of an n×m matrix; The study area is divided into k×k sub-areas; local statistical characteristics are calculated for each sub-area: mean μl=∑X(i, j) / n; standard deviation σl=√[∑(X(i, j)-μl) 2 / n]; calculate global statistical characteristics: mean μg=∑X(i,j) / N; standard deviation σg=√[∑(X(i,j)-μg) 2 / N]; Generate decomposition coefficient matrix.

[0086] Hydrological similarity calculation, including: defining characteristic vector F(i, j) = [H(i, j), S(i, j), A(i, j)]; where S is the slope and A is the catchment area; calculating similarity index SI(i, j) = exp(-||F(i, j) - F(i±1, j±1)|| / σ); generating similarity index; Identify outliers, including: setting a threshold θ = μ ± 2σ; marking the outlier position: if |X(i, j) - μ| > θ, then M(i, j) = 1, otherwise M(i, j) = 0; generate data quality indicators.

[0087] Filter data based on data quality indicators: If Q(i, j)>Qthreshold, retain the data. Otherwise, mark it as an invalid value; generate a valid data set.

[0088] Step S12: Specific process of hydrological process feature extraction Calculate the slope, S(i,j)=arctan[(H(i+1,j)-H(i-1,j)) 2 / 2dx+(H(i,j+1)-H(i,j-1)) 2 / 2dy]; Calculate the slope aspect A(i, j) = arctan[(H(i, j+1)-H(i, j-1)) / (H(i+1, j)-H(i-1, j))]; Obtaining geomorphic unit features and generating geomorphic feature vectors; Basis function construction, including: Construction of spatial distribution function: φ(x, y) = exp(-[(xx o ) 2 +(yy o )2 ] / 2σ 2 );Calculate the confluence path length: L=∑√[(xi+1-xi) 2 +(yi+1-yi) 2 ]; Generate basis function set B; Time-delay response analysis, including: calculation of the mutual correlation coefficient R(τ)=∑[P(t)*Q(t+τ)] / [σp*σq]; Identify the maximum correlation lag τmax=argmax(R(τ)); generate the lag coefficient matrix; Feature projection: Construct projection matrix: W=B*B T ; Calculate the transformation coefficient: T C =W*V; generate feature component set; Multi-scale combination: set the scale sequence s=[s1, s2, ..., sn]; calculate the energy of each scale E(s)=∑|F(s)|2; generate the hydrological feature set and process feature set; The data quality index is an n×m matrix with a value range of [0, 1]; the valid data set is a multidimensional data set that has passed the quality inspection; the landform feature vector is a triplet of [slope, aspect, elevation] for each grid point; the basis function set B is a set of orthogonal basis functions {φ i (x, y)}; the time lag coefficient matrix is ​​a t×t mutual correlation coefficient matrix; the characteristic component set F is a k-dimensional characteristic vector set; the hydrological characteristic set and the process characteristic set: a multidimensional characteristic matrix.

[0089] Step S2: Dynamic zoning and scale identification of mountain watersheds Hydrological feature set HF: a three-dimensional matrix of n×m×p, containing p hydrological feature layers; process feature set PF: a t×q matrix, containing q process features at t time points; Step S21: Specific process of dynamic partitioning algorithm Calculate terrain gradient, including: calculate x-direction gradient Gx(i, j) = [H(i+1, j) - H(i-1, j)] / 2Δx; calculate y-direction gradient Gy(i, j) = [H(i, j+1) - H(i, j-1)] / 2Δy; synthesize total gradient G(i, j) = √ [Gx(i, j) 2 +Gy(i,j) 2 ]; Generate gradient matrix G; Construct spatial neighborhood relationship: define the search radius r; calculate the spatial distance matrix D(i, j, k, l) = √[(xi-xk) 2 +(yj-yl) 2]; construct the neighborhood connection strength C(i, j, k, l) = exp(-D(i, j, k, l) / r)•G(i, j) / G(k, l); generate the neighborhood matrix N; Confluence time calculation: Calculate flow velocity V(i, j) = α•S(i, j) β •A(i, j) γ ; Calculate the time to the exit T(i, j) = ∑[L(i, j) / V(i, j)]; Generate the time matrix T; Step S22: Specific process of multi-scale decomposition Characteristic function construction: Construct terrain modulation function M(x, y) = exp(-κ•G(x, y)); define scale function: ψ(x, y, s) = M(x, y)•exp[-(x 2 +y 2 ) / 2s 2 ]; Generate characteristic function set F; 2. Multidimensional decomposition: Set the scale sequence S = {s1, s2, ..., s k}; Calculate the decomposition coefficient C (i,j,k)= JJf(x,y)·ψ(x,y,s k )dxdy; generate the decomposition coefficient set S; JJ is the double integral.

[0090] Step S3: Ecological flow pattern identification and uncertainty analysis Read the multi-scale feature matrix M, an n×m×k matrix containing features of k scale layers; Step S31: Specific process of eco-hydrological sequence matching Extraction of ecological response characteristics: calculate the flow change rate ΔQ(t)=[Q(t)-Q(t-1)] / Q(t-1); calculate the ecological sensitivity ES(t)=ΔQ(t) / ΔE(t) where ΔE(t) is the rate of change of environmental variables; generate the response feature set R.

[0091] The water system connectivity analysis includes: constructing a river section connection matrix cN(i, j) = {1, if river sections i and j are directly connected {0, otherwise; calculating the connectivity index cI(i) = ∑[cN(i, j) • W(j)] where W(j) is the weight of river section j; generating a connectivity matrix; Similarity calculation, including: constructing feature distance function d(x, y)=√[∑(w i •(x i -y i ) 2 )] where w i is the feature weight; calculate the time series similarity S(i, j)=exp[-d(R i , R j) / σ]; Generate similarity index S; Step S32: Specific process of uncertainty transfer analysis; Transfer network construction: define the node set V = {v i}, i=1,2,...,n; construct the adjacency matrix A(i,j)={w(i,j), if nodes i,j are connected {0, otherwise; generate the network structure N; Sensitivity analysis, including: calculation of local sensitivity LSI(i)=dY / dX i ; Calculate the overall sensitivity TSI(i)=Var(E[Y|X i ]) / Var(Y); Generate sensitivity matrix S; d is the sign of partial derivative; Step S4: Intelligent estimation of ecological flow Step S41: Adaptive dynamic integration specific process Model weight calculation, including: calculation of model error E(i, t) = [Q * i (t)-Q(t)] 2 ; Calculate the dynamic weight W(i, t) = exp[-E(i, t) / λ] / ∑exp[-E(j, t) / λ]; Generate the weight vector W; Reliability assessment, including: calculating the prediction interval PI(t)=[Q*(t)±z•σ(t)]; calculating the confidence level CR(t)=P[Q(t)∈PI(t)]; generating the reliability matrix; Step S42: Specific process of multi-constraint optimization Constraint construction, including: minimum ecological demand constraint Q(t)≥Q e (t); Maximum available water quantity constraint Q(t)≤Qm(t); Generate constraint condition set.

[0092] Objective function construction, including: ecological objective f1(Q)=∑[w1*(Q(t)-Q e (t)) 2 ]; water use target f2(Q)=∑[w2*(Q(t)-Q u (t)) 2 ]; Generate objective function F; Case 2: In some embodiments, the calculation process of the sub-steps may also be: Step S121: Geomorphic unit feature extraction: read a valid data set V: a multidimensional data set containing n×m grids; Slope calculation, including: x-direction slope Sx(i, j)=[H(i+1, j)-H(i-1, j)] / 2Δx; y-direction slope Sy(i, j)=[H(i, j+1)-H(i, j-1)] / 2Δy; composite slope: S(i, j)=√[Sx(i, j)2+Sy(i, j)2]; generate slope matrix S1; Slope calculation, including: slope angle A(i, j) = arctan[Sy(i, j) / Sx(i, j)]; direction adjustment A c (i, j) = {A (i, j), Sx (i, j) > 0 {A (i, j) + π, Sx (i, j) < 0; generate the aspect matrix S2; Elevation extraction, including: neighborhood average, Hm(i, j)=∑H(i±1, j±1) / 9; anomaly test, if|H(i, j)-Hm(i, j)|>θh: H(i, j)=Hm(i, j); generate elevation matrix S3; Δx, Δy grid spacing, A(i, j) original slope angle at position (i, j), A c (i, j) is the corrected slope angle at position (i, j), radian Hm(i, j) is the average elevation of the neighborhood at position (i, j), θh is the elevation anomaly detection threshold.

[0093] Step S122: constructing a spatial feature function, reading the topographic feature vector G as an n×m×3 matrix, including three feature layers: slope, aspect, and elevation.

[0094] The process of constructing the spatial distribution function: establish a local coordinate system based on the position (i, j) and calculate the spatial distribution function SP(x, y) = exp[-(x-xi) 2 / 2σx 2 -(y-yi) 2 / 2σy 2 ]·G(i, j), where σx and σy are characteristic diffusion parameters, and a set of spatial functions is generated. The D8 algorithm is used to track the flow path from position (i, j) to the outlet, and the path length is calculated as L(i, j)=∑√[(xk+1-xk)2+(yk+1-yk)2], where (xk, yk) are the coordinates of the sequence points on the path, and the path length matrix D is generated. Combining the spatial distribution function and the path length, the basis function B(x, y)=SP(x, y)·exp[-L(x, y) / L o ], where L ois the characteristic length parameter, and a basis function set is generated. Extract the rainfall sequence PR(t) and runoff sequence QR(t) from the valid data set, t=1, 2, ..., T, and generate the rainfall sequence and runoff sequence respectively. Calculate the rainfall-runoff correlation coefficient C(τ)=∑[PR(t)·QR(t+τ)] / [σp·σq], where τ is the time lag, and generate the correlation coefficient matrix. Identify the time lag τmax=argmax(C(τ)) corresponding to the maximum correlation coefficient, and calculate the response time of each sub-basin to generate the time lag coefficient matrix.

[0095] SP(x, y) is the spatial distribution function value at position (x, y), σx, σy are the characteristic diffusion parameters of the spatial distribution function, L(i, j) is the length of the confluence path from position (i, j) to the outlet, and L o is the characteristic length parameter of the basis function, PR(t) is the rainfall intensity at time t, QR(t) is the runoff at time t, and C(τ) is the correlation coefficient under time lag τ.

[0096] Step S123: Feature transformation and adjustment Basis function set B: p×q matrix, containing q basis functions; time delay coefficient matrix R: n×m time delay distribution matrix; effective data set V: n×m×k multidimensional data matrix; according to the orthogonalization of the basis function, the projection operator W(i, j)=∑[B(i, k)·B(j, k)] / λk is obtained, where λk is the eigenvalue of the kth basis function, and the projection matrix W is generated. Perform transformation operation T on the effective data C (i, j) = ∑ [W (i, k) · V (k, j)], where i = 1, 2, ..., n, j = 1, 2, ..., m, generates a transform coefficient set T C . Adjust according to the time lag coefficient A C (i, j) = T C (i, j)*exp[-R(i, j) / τ o ], where τ o is the characteristic time parameter, and generates an adjustment coefficient set. The adjustment coefficient is decomposed into F(i)=∑[A C (i, j)*ek(j)], where ek is the eigenvector, and the first k principal components are taken to generate the eigencomponent set.

[0097] Step S124: Multi-scale decomposition and reorganization The feature component set is an n×k matrix containing k principal components; set the multi-scale window length sequence K(i)={2 i |i=1, 2, ..., m}, where m is the maximum number of decomposition layers, and a window sequence K is generated. Wavelet decomposition D is performed on each scale C (i,j,k)=∑[F(t)·ψ((t-iτ) / s j)] / √s j , where ψ is the wavelet basis function, s j is the scale parameter, and a set of decomposition coefficients is generated. Calculate the energy distribution of each scale E(j)=∑|D C (i, j, k)| 2 , after normalization, the energy distribution matrix is ​​generated. According to the energy cumulative contribution rate, the important scale TV(j)=min{k|∑E(i) / ∑E(t)>η, i=1,...,k} is determined, where η is the energy threshold, and the threshold vector is generated. The important scale is selected for reconstruction HF(t)=∑[D C (i, j, k)·ψ((t-iτ) / s j )], generate the hydrological feature set. Based on the reconstruction results, calculate the feature parameter PF(i)=Φ[HF(t)], where Φ is the feature extraction operator, and generate the process feature set.

[0098] In another embodiment of the present application, step S211: terrain gradient calculation The hydrological feature set HF is an n×m×p matrix, containing p hydrological feature layers; the process feature set PF is a t×q matrix, containing q process features; The second-order central difference format is used to calculate the x-direction gradient Gx(i, j)=[H(i+1, j)-H(i-1, j)] / 2Δx and the y-direction gradient Gy(i, j)=[H(i, j+1)-H(i, j-1)] / 2Δy, where H is the elevation value, Δx and Δy are the grid spacing, and the height difference matrix H is generated. The composite slope S(i, j)=arctan(√[Gx(i, j)2+Gy(i, j)2]) is calculated based on the directional gradient to generate the slope matrix. Combine the elevation gradient and slope to calculate the composite gradient G(i, j)=√[Gx(i, j) 2 +Gy(i,j) 2 ]·[1+w·S(i,j)], where w is the slope weight coefficient, generates the gradient matrix G.

[0099] Step S212: Spatial neighborhood construction Read the gradient matrix, n×m gradient distribution matrix Set the search radius r, and identify each grid (i, j) that satisfies the distance condition d = √[(x-xi) 2 +(y-yj) 2 ]≤r, generating a neighborhood set N C . The connection strength K(i, j, k, l) = exp(-d / r o )·exp(-|G(i,j)-G(k,l)| / g o ), where r o is the distance characteristic parameter, go is the gradient feature parameter, and generates the connection strength matrix K. Set the threshold θ according to the connection strength c , when K(i, j, k, l)>θ c When the neighborhood connection relationship is established, the neighborhood matrix is ​​generated.

[0100] Step S213: Calculation of confluence time Read the neighborhood matrix: n×m×k×l neighborhood relationship matrix; read the gradient matrix: n×m gradient distribution matrix.

[0101] Based on the cumulative flow threshold θf, the river network unit is identified, the river network direction and connection relationship are extracted, the river network topology RN(i, j)={1: river network, 0: non-river network} is constructed, and the river network matrix RN is generated. The Manning formula is used to calculate the river flow velocity V(i, j)=(1 / n)·R(i, j) 2 / 3 S(i, j) 1 / 2 , where n is the roughness coefficient, R is the hydraulic radius, and S is the hydraulic slope, generating the velocity matrix V. Based on the D8 algorithm, the confluence path from each unit to the outlet is tracked, and the path length and cumulative time P(i, j)={(x1, y1, t1),...,(x n ,y n , t n )}, generate the path matrix P. Accumulate the transmission time of each path segment T(i, j) = ∑[L(k) / V(k)], where L(k) is the length of the kth path, and generate the time matrix.

[0102] Step S214: Dynamic partitioning rule construction Read the neighborhood matrix: neighborhood relationship matrix; read the time matrix T: confluence time distribution matrix.

[0103] Design space constraint function F(i, j) = α·N(i, j)+β·exp[-|T(i)-T(j)| / τ c ], where α, β are weight coefficients, τ c Generate constraint function for time characteristic parameter.

[0104] Based on the terrain characteristics and hydrological responses, the similarity between units is calculated as SM(i, j) = F(i, j) exp[-∑wk(xki-xkj)2], where wk is the feature weight and xk is the eigenvalue, and the similarity matrix SM is generated.

[0105] The hierarchical clustering method is used to determine the number of partitions and the similarity threshold TH={θ1, θ2, ..., θn}, and generate the threshold vector TH.

[0106] The partitioning criterion c(i, j) = {1: SM(i, j) > θk, 0: other} is constructed by combining spatial constraints and similarity thresholds to generate partitioning criteria. Among them, α and β are the weight coefficients of the constraint function, τ c is the time feature parameter, unit: hour, wk is the weight of the kth feature, θk is the similarity threshold of the kth partition, Step S215: Partition optimization and correction Partition identification initialization process: randomly assign an initial identification ID(i, j)=k, k∈[1, K] to each unit, where K is the preset number of partitions, and generate an identification matrix. Update the identification based on the objective function J=∑[w1·D(i)+w2·c(i, j)], where D(i) is the distance within the group, w1 and w2 are weight coefficients, and get oM(i, j)=argmin(J), and generate the optimization matrix oM. Calculate the spatial continuity index of each partition C I(k)=N c (k) / N(k), where N c (k) is the number of connected units, N(k) is the total number of units, and the continuity index cI is generated. The areas that do not meet the continuity requirements are merged or split, I(i, j) = f[oM (i, j), cI(k)], and the initial partition result is generated.

[0107] Step S216: Seasonal dynamic adjustment Read the initial partition results, calculate the hydrological characteristics of each partition in different seasons SF(k, s) = ∑[w(i) · X(i, s)], where w(i) is the feature weight and s is the seasonal index, and generate the seasonal feature matrix SF.

[0108] Design an adjustment function based on seasonal changes AF(t)=a o +∑[a i ·cos(2πit / T)+b i ·sin(2πit / T)], where T is the period, to generate the adjustment function AF. The dynamic range of the partition boundary is calculated based on the adjustment function CR(k, t)=r o ·[1+λ·AF(t)], where r o is the reference range, λ is the adjustment coefficient, and the change range matrix CR is generated. Combined with the initial partition and boundary changes, dynamic adjustment D(i, j, t) = f[I(i, j), cR(k, t)] is performed to generate a dynamic response unit.

[0109] Step S221: constructing a characteristic function and reading a dynamic response unit.

[0110] Terrain feature extraction process: calculate the morphological features of each unit TF(i) = {H(i), S(i), A(i), TWI(i)}, where H is elevation, S is slope, A is area, and TWI is terrain wetness index, and generate a terrain feature set TF.

[0111] Based on rainfall runoff data, the response characteristic RP(i)={t c (i), K c (i), η(i)}, where t c is the confluence time, K c is the flow generation coefficient, η is the flow generation efficiency, and the response parameter set RP is generated.

[0112] The characteristic expression EM(i, j) = φ(TF(i))·ψ(RP(i)) is constructed through orthogonalization, where φ and ψ are characteristic mapping functions, and an expression matrix is ​​generated.

[0113] Construct the characteristic function based on the expression F(x, y) = ∑[EM(i, j) · B(x, y)], where B is the basis function, and generate the characteristic function set F.

[0114] Step S222: Multidimensional decomposition implementation Read the characteristic function set. Design a multidimensional orthogonal basis function BM(x, y, z) = ∑[αk·Pk(x)·Qk(y)·Rk(z)], where P, Q, R are one-dimensional orthogonal polynomials, and generate the basis function matrix BM.

[0115] Construct the geometric scale sequence SS(k)=s o ·q k , k = 0, 1, ..., N, where s o is the initial scale, q is the scale ratio, and a scale sequence is generated. Multi-scale projection is performed on the input function ID(i, j, k) = ∬ F(x, y) · BM(x, y, k) dxdy to generate an initial decomposition set. Sparse optimization is performed using the L1 norm constraint S(i, j, k) = argmin{||F-BM·S|| 2 2+μ||S||1}, generating the decomposition coefficient set S.

[0116] Step S223: Energy distribution calculation Read the decomposition coefficient set S: Calculate the energy EV(k)=∑|S(i, j, k)| for the decomposition coefficients of each scale. 2 / N, where N is the total number of samples, to generate an energy value matrix. The kernel density estimation method is used to construct the energy distribution function DF(k)=(1 / nh)∑K[(x-xi) / h], where h is the bandwidth parameter and K is the kernel function to generate the distribution function. The relative energy contribution E(k)=EV(k) / ∑EV(i) of each scale is calculated and the cumulative energy analysis is performed to generate an energy matrix.

[0117] Step S224: Dominant mode recognition Read the energy matrix E, set the energy threshold and spatial correlation threshold IR={θe, θr}, where θe is the energy contribution rate threshold and θr is the correlation coefficient threshold, and generate the recognition criterion IR. Based on the recognition criterion, screen the dominant mode FM(k)={M(i)|E(i)>θe, R(i)>θr}, where R(i) is the spatial correlation coefficient, and generate the feature matrix FM. Calculate the spatiotemporal distribution characteristics of the dominant mode DM(i, j, k)=φ[FM(k)]·ψ(i, j), where φ and ψ are mapping functions, and generate the distribution matrix DM. Integrate the features and distribution information P(i, j, k)=γ·FM(k)+(1-γ)·DM(i, j, k), where γ is the weight coefficient, and generate the pattern matrix P.

[0118] Step S225: Scale information flow analysis Read the pattern matrix P, design the scale transfer function TF(s1, s2)=∫ψ(s1t)·ψ*(s2t)dt based on wavelet transform, where ψ is the wavelet function, and generate the transfer function TF. Calculate the information transfer amount TM(i, j)=∑P(i, k)·TF(sk, sj)·P(j, k), where sk and sj are scale parameters, and generate the transfer matrix TM. Evaluate the information flow intensity SI(k)=-∑TM(i, k)·l based on information entropy o g[TM(i, k)], generates the intensity index SI. The comprehensive transfer matrix and intensity index L(i, j) = TM(i, j) · [SI(i) + SI(j)] / 2, generates the flow index L.

[0119] Step S226: Feature Recombination Implementation Read the pattern matrix P and the flow index L. Calculate the reorganization weight W(k)=L(k)·E(k) / ∑[L(i)·E(i)] based on the information flow and energy distribution, and generate the weight vector W. Perform the feature reorganization operation RM(i, j)=∑[W(k)·P(i, j, k)], where k is the pattern number, and generate the reorganization matrix RM. Calculate the information retention rate RI(k)=||RM(k)-P(k)|| / ||P(k)|| of the reorganization result, and generate the representative index RI. Adjust the weight according to the representative index and recalculate M(i, j)=f[RM(i, j), RI(k)] to generate the multi-scale feature matrix M.

[0120] Step S311: Ecological response feature extraction Read the multi-scale feature matrix M, calculate the flow change rate of adjacent time periods WS(t)=[Q(t)-Q(t-1)] / Q(t-1), where Q(t) is the flow at time t, and generate the water volume sequence WS. Based on the characteristics of ecological water demand, construct the sensitivity function SP(t)=α·B(t)+β·E(t), where B(t) is the basic water demand and E(t) is the ecological water demand, and generate the sensitive period matrix SP. Calculate the ecological water demand for each period DM(t)=max[Qmin(t,Qe co (t)], where Qmin is the minimum water requirement, Qe co The demand matrix DM is generated for ecological water demand. The intensity matrix IM is generated by integrating water quantity changes and ecological demand IM(t)=WS(t)·SP(t) / DM(t). The response feature set R is generated by integrating various indicators R(t)={WS(t), SP(t), DM(t), IM(t)}.

[0121] Step S312: Water system connectivity analysis Read the multi-scale feature matrix M. Use the threshold method to extract the river network structure RN(i, j) = H[F(i, j) - θf], where H is the step function and θf is the extraction threshold, to generate the river network matrix RN. Calculate the connectivity P based on the topological relationship of the river section C (i,j)=d(i,j)·exp[-L(i,j) / L o ], where d is the connectivity indicator function, L is the length of the river section, and the physical connection matrix P is generated. C . Analyze the flow transmission characteristics H C (i, j) = Q (i, j) / [Q (i) + Q (j)], where Q is the cross-sectional flow, generating the hydrological connection matrix H C . Evaluate the ecological corridor function of river sections C I(i)=w1·P C (i)+w2·H C (i) Generate corridor index matrix CI. Integrate indicators C (i,j)=f[P C (i, j), H C (i, j), CI(i)], generating the connectivity matrix C .

[0122] Step S313: Similarity metric construction Read the response feature set R, connectivity matrix C , design the feature weight function WF(i)=exp[-λ·V(i)], where V(i) is the feature variation coefficient and λ is the attenuation parameter, and generate the weight function WF. Calculate the correlation coefficient between features C M(i, j) = cov(Xi, Xj) / [σ(Xi)·σ(Xj)], where cov is the covariance and σ is the standard deviation, generating the correlation matrix CM. Based on the Mahalanobis distance, the distance function DF(x, y) = √[(xy) T Σ -1 (xy)], where Σ is the covariance matrix, and the distance function is generated. The similarity S(i, j) = WF(i) exp[-DF(xi, xj) / h] is calculated by combining the weight and distance function, where h is the bandwidth parameter, and the similarity index is generated.

[0123] Step S314: Time series decomposition Read the multi-scale feature matrix and set the time window sequence WS(k)={w1,w2,...,w k}, where w i The window length is generated, and the window sequence is generated. Construct the decomposition function DF(t, s)=∫x(τ)ψ[(t-τ) / s]dτ based on wavelet transform, where ψ is the wavelet basis function, and generate the decomposition function.

[0124] Perform multi-scale decomposition DR(t, s) = DF(t, s) M(t), where s is the scale parameter, to generate the decomposition result.

[0125] The reorganization parameters are optimized according to the energy distribution, and T(t)=∑[α(s)·DR(t,s)] is obtained to generate a time series set.

[0126] Step S315: Pattern matching implementation Read similarity index S and time series set T. Design dynamic matching criteria M C (i, j) = S (i, j) · exp [-|t (i) - t (j)| / τ], where τ is the time characteristic parameter, generates the matching criterion.

[0127] The dynamic programming algorithm is used to calculate the sequence similarity DS(i, j) = min{d(i, j) + min[DS(i-1, j), DS(i, j-1), DS(i-1, j-1)]}, where d is the local distance, to generate a dynamic similarity matrix.

[0128] Based on the optimal path, feature matching points MP(k)={(i, j)|DS(i, j)<θd, path(i, j)=1} are extracted, where θd is the distance threshold, and a matching point set is generated. The matching results are integrated to generate a traffic pattern P(t)=f[MP(k), DS(i, j)], and a traffic pattern set is generated.

[0129] Step S321: Transmit network construction Read the traffic pattern set P, identify the key process nodes NP(i)={n|P(n)>θp}, where θp is the node threshold, and generate a node set. Construct a connection matrix CM(i, j)=G(Xi→Xj) based on causal tests, where G is the Granger causal test statistic, and generate a connection matrix. Design the traffic transfer function TF(i, j)=α·exp[-β·t(i, j)], where t is the transfer time, and generate a transfer function. Integrate node and connection information N(i, j)=CM(i, j)·TF(i, j) to generate a network structure.

[0130] Step S322: Sensitivity analysis implementation Read the traffic pattern set P and network structure N. Construct the parameter perturbation sequence DP(i)=P o (i)[1+δ·ξ(i)], where δ is the disturbance intensity and ξ is the random disturbance, and the disturbance scheme DP is generated. For each set of disturbances, the response result TR(i, j)=M[DP(i)] is calculated, where M is the model response function, and the test result TR is generated. The Sobol method is used to calculate the sensitivity index SM(i)=Var[E(Y|Xi)] / Var(Y), where Y is the response variable, and the sensitivity matrix SM is generated. Based on the sensitivity ranking, the importance score IV(i)=w1·SM(i)+w2·R(i) is calculated, where R is the correlation index, and the importance vector IV is generated. The sensitivity feature S(i, j)=f[SM(i), IV(j)] is generated by integrating various indicators, and the sensitivity matrix S is generated.

[0131] Step S323: Transfer path analysis Read the network structure N and sensitivity matrix S. Design the path identification criterion IR(p)=∏S(i)·∏N(i,j), where p is the path sequence, and generate the identification criterion IR. Calculate the path propagation characteristics PM(i,j)=∑[w(k)·φ(x(k))], where φ is the characteristic function, and generate the propagation characteristic matrix. Screen the critical path KP(i)={p|IR(p)>θr, PM(p)>θp} based on the feature importance, and generate the critical path set. Calculate the path importance index IM(i)=KP(i)·S(i) / ∑S(j), and generate the importance matrix. Integrate the path information H(i)=f[KP(i),IM(i)] and generate the transmission path set.

[0132] Step S324: Coupling effect evaluation Read the transfer path set H, determine the uncertainty source SM(i, j)=∑H(k)·U(k, i, j) based on path tracing, where U is the uncertainty distribution function, and generate the source matrix. Design the coupling effect evaluation function EF(i, j)=α·SM(i, j)+β·I(i, j), where I is the interaction intensity, and generate the evaluation function EF. Analyze the uncertainty interaction IM(i, j)=EF(i, j)·[1+γ· C (i, j)], where C is the coupling coefficient, generating the interaction matrix IM.

[0133] Action mode recognition process: Based on the interaction features, construct the pattern PM(i, j) = ψ[IM(i, j)], where ψ is the pattern recognition function, and generate the action mode matrix. Integrate the interaction features C (i)=f[IM(i),PM(i)], generating a set of coupling coefficients.

[0134] Step S325: Comprehensive uncertainty assessment Read the transfer path set and coupling coefficient set, design the comprehensive evaluation function EF(x)=∑[w(i)·x(i)]·c(i), where w is the weight coefficient, and generate the evaluation function. Analyze the contribution of each path C M(i, j)=H(i) C (j) / ∑[H(k)· c(k)], generate the contribution matrix. Calculate the overall uncertainty level LI(t)=∑[cM(i, t)·EF(i)], and generate the level index LI. Construct the spatiotemporal distribution feature DM(i, j)=φ[LI(i)]·ψ(j), where φ, ψ are mapping functions, and generate the distribution matrix. Comprehensively evaluate the results U(i)=f[LI(i), DM(i)], and generate the uncertainty index set.

[0135] Step S411: Extracting watershed characteristic parameters Read the flow pattern set and uncertainty index set. Extract the DEM derived parameters TF(i, j) = {H(i, j), S(i, j), A(i, j), TWI(i, j)}, where TWI is the terrain wetness index, and generate the terrain feature matrix.

[0136] Calculate the hydrological response characteristics RI(i, t) = Q(i, t) / [P(t) · A(i)], where Q is flow, P is rainfall, and A is area, and generate the response index matrix RI. Statistical land use composition LF(i) = {f1(i), f2(i), ..., f n (i)}, where f k The land characteristic matrix is ​​generated by taking the proportion of the kth type of land as the land characteristic matrix. The climate variability index cI(t)=[P(t)-P*] / σp is calculated, where P* is the average rainfall and σp is the standard deviation, and the climate indicator matrix is ​​generated. The characteristic indicators F(i)=g[TF(i), RI(i), LF(i), cI(i)] are integrated to generate the characteristic indicator set.

[0137] Step S412: Model selection criteria construction Read the feature index set F, design the model evaluation function EF(m)=∑[w(i)·e(i,m)], where e is the evaluation index and w is the weight, and generate the evaluation function. Calculate the index weight IW(i)=[1-H(i)] / ∑[1-H(j)] based on the entropy weight method, where H is the information entropy, and generate the index weight matrix. Use the clustering method to determine the selection threshold TV(k)=μ k +σ k , where μ k , σ k Generate a threshold vector for the inter-class statistical parameter. Generate a selection rule by integrating the evaluation function and the threshold SR(m)={1: EF(m)>TV(k), 0: others}. Generate a selection criterion by integrating the rule and the weight K(m)=f[SR(m), IW(m)].

[0138] Step S413: Spatiotemporal heterogeneity analysis Read the feature index set, divide the unit cU(i, j)=h[F(i, j)] based on the feature similarity, where h is the classification function, and generate the calculation unit set cU. Use Moran's I to calculate the spatial autocorrelation SH(i, j)=[n·∑∑wij(xi-x*)(xj-x*)] / [S o ∑(xi-x*) 2 ], generate the spatial heterogeneity matrix. Calculate the time series variation index TM(i, t) = C V[x(i, t)], where cV is the coefficient of variation, generates a time variation matrix. Analyze multi-scale features SI(i, s) = x(i, s) / x(i, so ), where s is the scale parameter, and the scale indicator matrix is ​​generated. The spatiotemporal features H(i, j) = w1*SH(i, j) + w2*TM(i, j) + w3*SI(i, j) are integrated to generate the heterogeneity matrix.

[0139] Step S414: Weight allocation optimization reads selection criteria and heterogeneity matrix. Design weight allocation function DF(x)=exp[-λ·x] / ∑exp[-λ·x], where λ is the temperature parameter, and generate allocation function DF. Calculate model adaptability score SM(i, j)=K(i)·[1-H(j)], where K is the selection score, and generate score matrix. Construct optimization objective function O P(w)=max{∑[SM(i)·w(i)]|∑w(i)=1}, generate the optimization solution. Solve the optimization solution to get the weight vector W(i)=DF[OP(i)], generate the weight vector.

[0140] Step S415: Reliability assessment is implemented by reading the traffic pattern set, constructing a cross-validation scheme TP(i)={T(i), V(i)}, where T is the training set and V is the validation set, and generating a test scheme. Perform model validation VR(i, j)=M[T(i), V(j)], where M is the model operation, and generate a validation result VR. Calculate multiple precision indicators AM(i)={NSE(i), RMSE(i), R2(i)} to generate a precision matrix. Evaluate the model stability index SI(i)=σ[VR(i)] / μ[VR(i)] to generate a stability index. Comprehensively evaluate the index R(i, j)=f[AM(i), SI(j)] to generate a reliability matrix. Among them, NSE is the Nash efficiency coefficient, RMSE is the root mean square error; R2 is the determination coefficient, σ is the standard deviation operator, and μ is the mean operator.

[0141] Step S416: The integrated model constructs the read weight vector and reliability matrix. Construct a hierarchical integrated structure SF(i)={L1(i), L2(i), ..., L n (i)}, where L k Generate an integrated framework for the k-th model group. Design model combination function C S(m)=α·W(m)·R(m)+β·D(m), where D is the diversity index, to generate a combination strategy. Use cross-validation to optimize the integrated parameters PM(i, j)=argmin{∑[yt-f(xt, θ)]2} to generate a parameter matrix. Based on the optimal parameters, build an integrated model E(x)=∑[CS(i)·Mi(x)] to generate an integrated model E.

[0142] Step S421: Constraints are used to construct a read integration model E, and analyze the ecosystem requirements NF(t)={Q1(t), Q2(t), ..., Q n (t)}, where Q i Generate the demand characteristic matrix for the i-th type of demand flow. Set the upper and lower limit constraints W C (t)=[Qmin(t), Qmax(t)], generate water constraint matrix. Protection target setting process: Construct ecological protection index PM(i)=g[B(i), E(i)], where B is the reference value and E is the target value, generate protection target matrix. Calculate system carrying capacity CI(t)=h[S(t), R(t)], where S is supply and R is demand, generate carrying capacity index. Integrate various constraints C(i, t)=f[W C (t), PM(i), CI(t)], generating a set of constraints C .

[0143] Step S422: The objective function constructs the read constraint condition set C. Construct the periodic function SF(t)=ao+∑[a i cos(ω i t)+b i sin(ω i t)], where ω i is the frequency parameter, and generates the variation function SF. Design multi-objective optimization criteria OC (x)={f1(x),f2(x),...,f n (x)}, where f i Generate optimization criteria for the i-th objective function OC . Design target trade-off function TF(x)=∑[w i ·f i (x)], where w i is the weight coefficient, generating the trade-off function TF. Integrating seasonal changes and the trade-off function F(x, t) = TF(x) · SF(t) generates the objective function F. Where a o , a i , b i is the periodic function coefficient, ω i is the frequency parameter, unit: radian / day; f i is the ith sub-objective function.

[0144] Step S423: Optimize the solution to read the constraint condition set and the objective function. Construct the optimization algorithm parameter set PS(i)={p1(i), p2(i), ..., p n (i)}, where p i For the algorithm parameters, generate parameter sets. Perform optimization iterations IR(k)=argmin{F(x)| C(x)≤0}, where k is the number of iterations, and the iterative results are generated. Evaluate the optimization convergence characteristics C M(k)=|IR(k)-IR(k-1)| / |IR(k-1)|, generate the convergence index. Based on the convergence criterion, select the optimal solution I(x)=IR[k*], where k* is the optimal number of iterations, and generate the initial result.

[0145] Step S424: Ecological response assessment Read the initial results and design the ecological response evaluation function EF(x)=∑[α i ·r i (x)], where r i is the response index, α i is the weight coefficient, and the evaluation function is generated. Calculate the ecological index change VM(i, t) = [y(i, t) - y o (i, t)] / y o (i, t), where y o is the baseline value, and a change matrix is ​​generated. Analyze the response sensitivity SI(i)=dy(i) / dx(i)·[x(i) / y(i)] and generate the sensitivity index SI. Extract the response feature parameters FM(i, j)=g[VM(i), SI(j)] and generate the feature matrix FM. Integrate the evaluation results B(i)=h[SI(i), FM(i)] and generate the response index B. g, h are feature extraction functions. It should be noted that d is the symbol of partial derivative.

[0146] Step S425: Result correction implementation Read the initial results and response indicators. Design the correction function MF(x)=x+λ·δ(x), where λ is the correction coefficient and δ is the deviation function, and generate the correction function. Determine the adjustment parameter AM(i, j)=f[B(i), I(j)] based on the response indicator and generate the adjustment matrix AM. Perform the correction calculation MR(i)=MF[I(i)]·AM(i) to generate the correction result. Evaluate the correction effect EI(i)=|MR(i)-I(i)| / |I(i)| and generate the effect indicator. Determine the final result R(i)=k[MR(i), EI(i)] based on the effect indicator and generate the ecological flow estimation result.

[0147] Step S426: Reliability evaluation implementation Read the ecological flow estimation results and construct the evaluation index system IS(i)={q1(i),q2(i),...,q n (i)}, where q i Generate an index system for evaluation indicators. Calculate the reliability index RM(i, j)=p[R(i, j)≥R o (i, j)], where R oAs the target value, generate a reliability matrix. Evaluate the uncertainty of the result UI(i)=σ[R(i)] / μ[R(i)] and generate an uncertainty index. Analyze the adaptability characteristics AM(i, j)=f[RM(i), UI(j)] and generate an adaptability matrix. Comprehensively evaluate the result A(i)=g[RM(i), UI(i), AM(i)] and generate a reliability evaluation index.

[0148] The preferred embodiments of the present invention are described in detail above; however, the present invention is not limited to the specific details in the above embodiments. Within the technical concept of the present invention, various equivalent transformations can be made to the technical solutions of the present invention, and these equivalent transformations all belong to the protection scope of the present invention.

Claims

1. An intelligent estimation method of ecological flow in mountainous watersheds without data based on multi-scale fusion, characterized by: The following steps are involved: Collect basic data of the study area and preprocess it to obtain effective data sets; basic data include DEM data, remote sensing image data, meteorological reanalysis data and land use data; Extract hydrological feature sets and process feature sets based on valid data sets; Generate dynamic response units based on hydrological feature sets and process feature sets, and perform multi-dimensional hydrological process decomposition to obtain a multi-scale feature matrix; The multi-scale characteristic matrix is ​​matched with eco-hydrological sequences to obtain the flow pattern set, and the hierarchical uncertainty transfer analysis is performed to generate the uncertainty indicator set; An integrated model was constructed based on the flow pattern set and uncertainty indicator set, and multi-constraint optimization was performed to obtain the ecological flow estimation results.

2. The method according to claim 1, characterized in that Preprocessing includes: Obtain and calculate the terrain complexity index value of basic data to obtain the terrain decomposition factor; Based on the terrain decomposition factor, the basic data is locally decomposed and globally decomposed to generate a decomposition coefficient matrix; Calculate the hydrological similarity value according to the decomposition coefficient matrix to obtain the similarity index; Generate data quality indicators based on similarity indicators to obtain a valid data set.

3. The method according to claim 1, characterized in that The steps of extracting hydrological feature sets and process feature sets based on valid data sets include: Extract the slope value, aspect value and elevation value of the grid unit, calculate the geomorphic position index and obtain the position index matrix; The geomorphic unit types are divided based on the position index matrix, and the characteristic parameters of each type of unit are extracted to generate the geomorphic feature vector; Constructing eigendecomposition basis functions based on topographic feature vectors to obtain a basis function set; Calculate the time response relationship between rainfall data and runoff data to obtain the time lag coefficient matrix; According to the basis function set and the time-delay coefficient matrix, the effective data set is subjected to feature extraction to obtain the feature component set; The characteristic component sets are combined at multiple scales to generate hydrological characteristic sets and process characteristic sets.

4. The method according to claim 1, characterized in that The steps to generate a dynamic response unit include: Calculate the terrain gradient value of each point in the watershed to obtain the gradient matrix; Construct spatial neighborhood relationships based on the gradient matrix and generate a neighborhood matrix; Calculate the confluence time of each area and obtain the time matrix; A dynamic partitioning standard is constructed based on the neighborhood matrix and the time matrix to obtain the partitioning criterion; Perform adaptive partitioning according to the partitioning criteria to obtain the initial partitioning result; The initial partitioning results are dynamically adjusted according to seasonal changes to obtain dynamic response units.

5. The method according to claim 4, characterized in that The steps for multidimensional hydrological process decomposition include: Extracting the terrain features of each unit to obtain a terrain feature set; Calculate the hydrological response parameters of each unit to obtain a response parameter set; Construct the characteristic expression of unit scale to obtain the characteristic function set; Perform multi-dimensional decomposition based on the characteristic function set to obtain an initial decomposition set; Optimize the decomposition coefficients of the initial decomposition set to generate a decomposition coefficient set; The energy value of each decomposition scale is calculated based on the decomposition coefficient set to obtain a multi-scale feature matrix.

6. The method according to claim 1, characterized in that The steps of obtaining the flow pattern set by matching the eco-hydrological sequence include: Extract the water volume change sequence to obtain the water volume sequence; Calculate the characteristics of ecological sensitive periods to obtain the sensitive period matrix; Analyze ecological water demand to obtain the demand matrix; The response intensity is calculated based on the water quantity series and the sensitive period matrix to obtain the intensity matrix; Combining the demand matrix and the intensity matrix to construct response characteristics and generate a response feature set; Pattern recognition is performed based on the response feature set to generate a traffic pattern set.

7. The method according to claim 5, characterized in that The steps for performing hierarchical uncertainty transfer analysis to generate uncertainty indicator sets include: Construct the hydrological process transmission chain network to obtain the network structure; Calculate the parameter sensitivity index of each node to generate a sensitivity matrix; Analyze the uncertainty transfer path according to the network structure and sensitivity matrix to obtain the transfer path set; Evaluate the coupling effect of multi-source uncertainties to obtain a set of coupling coefficients; The comprehensive uncertainty is calculated based on the transfer path set and the coupling coefficient set to generate an uncertainty index set.

8. The method according to claim 1, characterized in that The steps to build an ensemble model include: Extract terrain feature parameters to obtain terrain feature matrix; Calculate the hydrological response index to obtain the response index matrix; Analyze land use characteristics to obtain land characteristic matrix; Assess climate variability to obtain a climate indicator matrix; Model selection criteria were constructed based on terrain characteristics matrix, response indicator matrix, land characteristics matrix and climate indicator matrix; Calculate the spatiotemporal heterogeneity index of each region to obtain the heterogeneity matrix; Determine the model weights according to the model selection criteria and the heterogeneity matrix to generate a weight vector; Evaluate the reliability index of each model to obtain the reliability matrix; The integrated model is constructed using the weight vector and reliability matrix.

9. The method according to claim 1, characterized in that The steps of multi-constraint optimization to obtain ecological flow estimation results include: Construct the seasonal variation function to obtain the variation function; Design the ecological demand constraint system to obtain the constraint condition set; Construct an objective function that takes seasonal changes into account; Perform optimization and solution based on the constraint condition set and the objective function to obtain the initial result; Analyze the response characteristics of the ecosystem to the initial results to obtain response indicators; The initial results are revised based on the response indicators to generate ecological flow estimates.

10. The method according to claim 5, characterized in that Ecohydrological sequence matching also includes: Extract river network structure data to obtain river network matrix; Calculate the physical connectivity between river sections to obtain the physical connectivity matrix; Analyze the hydrological connectivity to obtain the hydrological connectivity matrix; Evaluate the ecological corridor function to obtain the corridor index matrix; The physical connectivity matrix, hydrological connectivity matrix and corridor index matrix were integrated to generate the connectivity matrix; The connectivity matrix, together with the response feature set, is used for pattern recognition to generate a set of traffic patterns.

11. The method according to claim 6, characterized in that The steps to analyze the uncertainty transfer path include: Design a perturbation experiment plan to obtain a perturbation plan; Perform parameter perturbation test to obtain test results; Calculate the response sensitivity to obtain the sensitivity matrix; Analyze the importance of parameters to obtain the importance vector; Identify the sources of uncertainty and obtain the source matrix; Analyze the coupling action mode to obtain the action mode matrix; Generate a set of coupling coefficients based on the action mode matrix.

Images