A method and system for constructing a rainfall runoff model

Abstract

This invention discloses a method and system for constructing a rainfall-runoff model. The method involves collecting and preprocessing relevant daily sequence data within a watershed; constructing a hydrological model component library; pre-assembling the hydrological model structure; defining a parameter space; performing a collaborative search using genetic programming; simulating daily runoff using the optimal model structure and parameters obtained through the search and calibration; evaluating the model's generalization performance; and outputting the optimal model structure type, optimized parameter vector, calibration set, validation set performance indicators, and runoff simulation sequence. This invention can automatically generate rainfall-runoff models with adaptive model structures and interpretable hydrophysical mechanisms, effectively improving the automatic discovery capability of model structures and the interpretability of hydrophysical mechanisms, while reducing manual workload while avoiding overfitting.

Description

Technical Field

[0001] This invention relates to the field of interdisciplinary technology of hydrological modeling and artificial intelligence, and in particular to a method and system for constructing a rainfall-runoff model. Background Technology

[0002] Conceptual hydrological models express hydrological processes such as interception, evapotranspiration, infiltration, runoff generation, and confluence through mathematical relationships, and are important tools for watershed hydrological simulation and flood forecasting. Currently, dozens of conceptual hydrological models have been developed internationally, which vary significantly in terms of structural complexity, number of parameters, and applicable watershed types.

[0003] Traditional hydrological modeling methods typically involve hydrologists manually selecting a specific model structure based on watershed characteristics and empirical knowledge, then calibrating the model parameters using trial and error, optimization algorithms, or empirical formulas. However, these methods rely heavily on expert experience for model structure selection, making them highly subjective and potentially overlooking more suitable models. Automatic parameter calibration usually assumes a fixed model structure, failing to achieve coordinated optimization of structure and parameters. Furthermore, complex models often perform well on calibration sets but degrade on validation sets, posing a risk of overfitting. Finally, for specific watersheds, there is a lack of systematic frameworks for model structure selection and parameter optimization, resulting in low modeling efficiency.

[0004] In recent years, the international hydrological community has proposed the concept of "model structure uncertainty," emphasizing that structure selection and parameter calibration should be considered simultaneously during model building. Open-source hydrological modeling frameworks such as SuperflexPy provide flexible model building tools, but manual design of model structures is still required.

[0005] Genetic programming, as a type of symbolic regression algorithm based on Darwinian evolution, possesses the core advantage of being able to search a broader function space through genetic operations, automatically discovering predictive models with explicit mathematical expressions. This process can automatically calibrate model parameters while searching for the optimal model structure, effectively achieving synergistic results in structure discovery and parameter optimization.

[0006] However, in the field of rainfall-runoff simulation, genetic programming has not yet fully realized its core potential in "model structure discovery." Existing applications mainly focus on short-term runoff forecasting, using genetic programming as a symbolic regression tool to construct functional rainfall-runoff relationships for short-term runoff forecasts from time-series data such as rainfall and flow. While this application demonstrates the nonlinear fitting capabilities of genetic programming, it is essentially still at the level of data-driven function fitting and has failed to leverage its advantages in co-optimization of model structure search and automatic model parameter calibration.

[0007] Therefore, it is necessary to apply the "model structure discovery" capability of genetic programming to the construction of rainfall-runoff models, and to invent a method for intelligent selection and automatic parameter calibration of hydrological models. By constructing a hydrological model component library and pre-assembling candidate model structures, the genetic programming algorithm is used to automatically achieve model structure selection and automatic parameter calibration, thereby obtaining an optimal model with good predictive performance, physical interpretability and good fit. Summary of the Invention

[0008] Purpose of the invention: This invention provides a method and system for constructing a rainfall-runoff model, which can automatically generate rainfall-runoff models with adaptive model structure and interpretable hydrophysical mechanisms, effectively improving the automatic discovery capability of model structure and the interpretability of hydrophysical mechanisms, and reducing the amount of manual work while avoiding overfitting.

[0009] Technical solution: The present invention provides a method for constructing a rainfall-runoff model, comprising the following steps:

[0010] Step 1: Collect and preprocess relevant daily sequence data within the watershed;

[0011] Step 2: Construct a hydrological model component library;

[0012] Step 3: Pre-assemble the hydrological model structure;

[0013] Step 4: Define the parameter space;

[0014] Step 5: Use genetic programming for cooperative search;

[0015] Step 6: Use the optimal model structure and parameters obtained through search and calibration to simulate daily runoff, evaluate the model's generalization performance, and output the optimal model structure type, optimized parameter vector, calibration set, validation set performance index, and runoff simulation sequence.

[0016] Furthermore, in step 1, the daily rainfall sequence of each rain gauge station in the basin, the daily outlet cross-sectional flow sequence of the basin, and the daily net radiation sequence, daily soil heat flux sequence, daily average temperature sequence, daily wind speed sequence at 2m height, and daily actual water vapor pressure sequence used for the calculation of potential evapotranspiration in the basin are collected. The obtained sequences are then preliminarily processed and divided into calibration set and validation set.

[0017] Furthermore, in step 2, a hydrological model component library is constructed, which includes eight basic hydrological components, including unsaturated reservoirs, power function reservoirs, linear reservoirs, interception filters, runoff-generating reservoirs, confluence modules, unit lines, and tank units. Each type of component implements independent water balance calculation and status update functions.

[0018] Furthermore, in step 3, the hydrological model structure is pre-assembled. Based on the hydrological component library, nine pre-assembled model structures are constructed, including HBV simple structure, HBV series structure, HBV parallel structure, GR4J structure, HYMOD structure, M4 structure, two-layer water tank structure, three-layer water tank structure and four-layer water tank structure. Each structure defines a specific component connection topology and data flow direction.

[0019] Furthermore, in step 4, a parameter space is defined, and a parameter vector and its physical constraint range are established for each structure. The parameters include common parameters, HBV model structural parameters, GR4J model structural parameters, HYMOD model structural parameters, and TANK model structural parameters.

[0020] Furthermore, in step 5, cooperative search is performed using genetic programming, specifically including the following steps:

[0021] Step 51, Genetic Programming Population Initialization: Randomly generate individuals within the population. Each individual contains a randomly selected model structure type and a parameter vector uniformly sampled in the corresponding parameter space.

[0022] Step 52, Fitness Assessment: For each individual in the population, hydrological simulation is performed using precipitation, potential evapotranspiration and measured runoff data within the set, and the individual fitness is calculated and assessed.

[0023] Step 53, Individual Selection: A tournament selection strategy is adopted, in which a specified number of individuals are randomly selected from the population each time for comparison, and the one with the highest fitness is selected as the parent.

[0024] Step 54, Parent Crossover: For the two selected parent individuals, perform structural crossover and parameter crossover with predefined crossover probabilities. Structural crossover swaps the model structure types of the two parents with predefined probabilities. For parameter crossover, for each parameter required by the offspring, randomly select and inherit from the parameters of the two parents at the corresponding position.

[0025] Step 55, Individual Mutation: For offspring individuals generated by crossover, Gaussian mutation is performed on each parameter with a predefined mutation probability. The mutation amount follows a normal distribution with a mean of 0 and a standard deviation of σ×(parameter upper bound - parameter lower bound), where σ is the mutation intensity coefficient. When the parameter exceeds the boundary after mutation, a mirror reflection strategy is used to bring it back to the effective range.

[0026] Step 56, Elite Retention: Sort the individuals in the current population from high to low fitness, retain a predefined number of elite individuals to directly enter the next generation population, and replace the remaining individuals with the offspring generated by the crossover and mutation operation.

[0027] Step 57, Early Termination Judgment: Record the optimal fitness of each generation. If the optimal fitness of the population does not improve for a certain number of consecutive generations, then terminate the evolution.

[0028] Step 58, Iterative Update: Repeat the above steps until the predefined maximum number of generations is reached or the early stopping condition is triggered.

[0029] Furthermore, step 6, evaluating the model's generalization performance, specifically includes the following steps:

[0030] Step 61: Use the rate set to perform hydrological simulation on the generated population, set a simulation warm-up period to make the model state tend to stabilize, and the simulation results during the warm-up period are not included in the subsequent index calculation, thereby reducing the impact of the uncertainty of the initial state;

[0031] Step 62: Calculate the NSE, KGE, RMSE and VE indices for each individual;

[0032] Step 63: Introduce a complexity penalty term CP based on the BIC concept to control model complexity and avoid overfitting. The specific formula is shown below:

[0033]

[0034] Where: λ is the regularization strength; n is the total time step of the calibration set data; k is the number of parameters of the generated model;

[0035] Step 64: Calculate the final model fitness function and use it to evaluate fitness.

[0036] Accordingly, a rainfall-runoff model construction system includes: a data acquisition module, a component construction module, a collaborative search module, and a model output module. The data acquisition module collects and preprocesses relevant daily sequence data within the watershed. The component construction module constructs a hydrological model component library, pre-assembles the hydrological model structure, and defines the parameter space. The collaborative search module uses genetic programming to perform collaborative search. The model output module uses the optimal model structure and model parameters obtained through search and calibration to simulate daily runoff, evaluate the model's generalization performance, and output the optimal model structure type, optimized parameter vector, calibration set, validation set performance indicators, and runoff simulation sequence.

[0037] Beneficial effects: Compared with the prior art, the present invention has the following significant advantages: By constraining the search space of genetic programming to the hydrological model component library and its pre-assembled structure, the present invention realizes the automated discovery of model structure and the integrated calibration of parameters, improves the probability of finding the optimal solution, and reduces the modeling threshold; by introducing a complexity penalty term, the present invention can automatically balance the model accuracy and model complexity during the model optimization process, effectively select models with simplified structure and strong generalization ability, and reduce the risk of model overfitting. Attached Figure Description

[0038] Figure 1This is a schematic diagram of the method flow of the present invention.

[0039] Figure 2 This is a diagram showing the optimal model structure and parameters obtained by searching and calibrating the Red Creek watershed based on the genetic programming algorithm of this invention.

[0040] Figure 3 The figure shows the rainfall-runoff simulation results of the optimal model searched and calibrated in the Red Creek watershed according to this invention.

[0041] Figure 4 This is an evolutionary curve of the search and calibration model based on the genetic programming algorithm in the Red Creek watershed according to the present invention. Detailed Implementation

[0042] like Figure 1 As shown, a method for constructing a rainfall-runoff model includes the following steps:

[0043] Step 1: Collection and preprocessing of relevant daily sequences within the watershed: Collect daily rainfall sequences from each rain gauge station within the watershed, daily outlet cross-sectional flow sequences from the watershed, as well as daily net radiation sequences, daily soil heat flux sequences, daily average temperature sequences, daily wind speed sequences at 2m altitude, and daily actual water vapor pressure sequences for calculating potential evapotranspiration in the watershed. Perform preliminary processing on the acquired sequences and divide them into calibration sets and validation sets.

[0044] Step 2, Hydrological Model Component Library Construction: Establish 8 basic hydrological components, including unsaturated reservoirs, power function reservoirs, linear reservoirs, interception filters, runoff-generating reservoirs, confluence modules, unit lines, and tank units. Each component implements independent water balance calculation and status update functions.

[0045] Step 3: Pre-assembly of hydrological model structure: Based on the hydrological component library, construct 9 pre-assembled model structures, including HBV simple structure, HBV series structure, HBV parallel structure, GR4J structure, HYMOD structure, M4 structure, two-layer water tank structure, three-layer water tank structure and four-layer water tank structure. Each structure defines a specific component connection topology and data flow direction.

[0046] Step 4, Parameter Space Definition: Based on the characteristics of each model structure, establish a parameter vector and its physical constraint range for each structure. The parameters include common parameters, HBV model structure parameters, GR4J model structure parameters, HYMOD model structure parameters, and TANK model structure parameters.

[0047] Step 5: Use genetic programming for cooperative search;

[0048] Step 6, Performance Evaluation: Daily runoff simulation is performed using the optimal model structure and parameters obtained through search and calibration to evaluate the model's generalization performance; Optimal Model Output: Output the optimal model structure type, optimized parameter vector, performance indices of calibration set and validation set, and runoff simulation sequence.

[0049] Step 1 provides the data input. Steps 2-5 involve dynamic model building and evaluation after the calibration set data input is provided in Step 1, searching for the optimal model. Step 6 uses the optimal model structure and parameters output in Step 5 to simulate daily runoff on the validation set and evaluate the model's generalization performance.

[0050] A method for intelligent optimization and automatic parameter calibration of hydrological models is presented, taking the Red Creek watershed in Utah, USA as an example. The following process is implemented based on the R language platform, and the method process is described in detail with reference to the attached figures.

[0051] Step 1: Collect daily observation data from October 1, 1980 to September 30, 2014 for the Red Creek watershed (corresponding USGS site number: 02479300, data sourced from the Camels dataset https: / / ral.ucar.edu / solutions / products / camels). This includes rainfall P (mm / day), potential evapotranspiration PET (mm / day), and watershed outlet cross-sectional flow Q_obs (mm / day). The data is divided into a calibration period and a validation period, comprising 75% and 25% of the total data, respectively. A 30-day warm-up period is set.

[0052] Step 2, Hydrological Model Component Library Construction: Establish 8 basic hydrological components, including unsaturated reservoirs, power function reservoirs, linear reservoirs, interception filters, runoff-generating reservoirs, confluence modules, unit lines, and tank units. Each component implements independent water balance calculation and status update functions.

[0053] Step 21: Construction of Unsaturated Reservoir Components: Based on the HBV model, the actual evapotranspiration and runoff are calculated by simulating the water balance of the soil layer. The actual evapotranspiration adopts a nonlinear relationship, while the runoff adopts a power function relationship. The specific formulas are as follows:

[0054]

[0055]

[0056]

[0057] Where P(t) is the precipitation at time t (mm / day); PET(t) is the potential evapotranspiration at time t (mm / day); AET(t) is the actual evapotranspiration at time t (mm / day); S(t) is the soil water storage at time t (mm); S max C represents the maximum water storage capacity of the soil (mm). e is the evapotranspiration correction factor; m is the evapotranspiration shape parameter; β is the flow generation index; R(t) is the flow generation rate at time t (mm / day); Δt is the time step (day);

[0058] Step 22: Construction of Power Function Reservoir Component and Linear Reservoir Component: Construct a component to simulate the groundwater or river regulation process, where the outflow and storage volume have a power function relationship. When the power exponent is 1, the power function reservoir degenerates into a linear reservoir. The specific formula is as follows:

[0059]

[0060]

[0061] Among them, R in Q(t) is the inflow rate (mm / day); Q(t) is the outflow rate (mm / day); k is the outflow coefficient (day). -1 S(t) represents the water storage capacity (mm); α represents the outflow index;

[0062] Step 23: Construction of the interception filter component: Based on the GR4J model, this component can simulate vegetation interception and direct evapotranspiration processes. It subtracts the evaporated portion from P and PET, outputting net precipitation and net evapotranspiration. The specific formulas are as follows:

[0063]

[0064]

[0065]

[0066] Where P(t) is the precipitation at time t; PET(t) is the evapotranspiration at time t; remove(t) is the amount removed at time t; P net (t) represents the net precipitation at time t; PET net (t) represents the net evaporation at time t;

[0067] Step 24: Construction of the runoff-generating reservoir component: Based on the GR4J model, this component reflects the dual impact of soil moisture content on runoff generation and evapotranspiration. A nonlinear runoff-generating relationship is used to simultaneously calculate evapotranspiration, runoff generation, and infiltration. The specific formulas are as follows:

[0068]

[0069]

[0070]

[0071]

[0072] Where x1 is the runoff storage capacity parameter; P s (t) represents the outflow rate (mm / day); P r (t) represents the infiltration rate (mm / day); n i α is the infiltration coefficient; β are shape parameters;

[0073] Step 25: Construction of the confluence module component: Based on the GR4J model, it adopts a power function outflow relationship and includes a river water exchange term, which can simulate the mutual recharge between surface water and groundwater. The specific formula is as follows:

[0074]

[0075]

[0076]

[0077] Where x2 is the groundwater exchange parameter; x3 is the reservoir capacity parameter; F(t) is the groundwater exchange flux (mm / day); Q(t) is the reservoir outflow (mm / day); and I(t) is the inflow (mm / day). , For shape parameters;

[0078] Step 26, Unitline Component Construction: This component expresses the time delay from precipitation occurrence to the formation of watershed outlet runoff. It simulates the lag effect of watershed runoff convergence using a time-shifting method. The specific formula is as follows:

[0079]

[0080]

[0081] Where Q(t) is the runoff at the outlet section at time t (mm / day); P(t-τ) is the inflow at time t (mm / day); and n is the total time step of the calibration set data. This is the lag time parameter;

[0082] Step 27, Water Tank Unit Component Construction: Based on the TANK model, this component serves as the basic unit of the TANK model and simulates layered soil structure through multi-level series connection. The specific formula is as follows:

[0083]

[0084]

[0085]

[0086]

[0087] Where a1 is the outflow coefficient of the side hole (day) -1 h1 represents the range of side hole elevation (mm); a2 represents the bottom hole outflow coefficient (day). -1 ); Inflow rate (mm / day); Q side (t) represents the flow rate of the side orifice (mm / day); Q bottom (t) represents the flow rate of the bottom hole (mm / day);

[0088] Step 3: Pre-assembly of hydrological model structure: Based on the hydrological model component library created in Step 2, pre-assemble the hydrological model structure, specifically including HBV simple structure, HBV series structure, HBV parallel structure, GR4J structure, HYMOD structure, M4 structure, two-layer water tank structure, three-layer water tank structure and four-layer water tank structure.

[0089] Step 31: Pre-assembly of HBV simple structure: unsaturated reservoir → power function reservoir, including two stages of runoff generation and runoff collection, with 4 parameters, specifically including maximum soil water storage capacity Smax, evapotranspiration coefficient Ce, runoff generation index beta and runoff discharge coefficient k;

[0090] Step 32: Pre-assembly of HBV series structure: unsaturated reservoir → power function reservoir 1 → power function reservoir 2. The two-stage confluence can better simulate the different response characteristics of fast flow and base flow. There are 7 parameters, including maximum soil water storage capacity Smax, evapotranspiration coefficient Ce, runoff index beta, power index alpha1 of the first reservoir, power index alpha2 of the second reservoir, drainage coefficient k1 of the first reservoir, and drainage coefficient k2 of the second reservoir.

[0091] Step 33: Pre-assembly of HBV parallel structure: Unsaturated reservoir → distributor → fast power function reservoir + slow power function reservoir, explicitly distinguishing between fast and slow runoff paths. There are 6 parameters, specifically including maximum soil water storage capacity Smax, evapotranspiration coefficient Ce, runoff index beta, fast path diversion ratio frac_fast, fast reservoir drainage coefficient k_fast, and slow reservoir drainage coefficient k_slow. The specific connection topology is as follows:

[0092] Rainfall and potential evapotranspiration are input into unsaturated reservoirs. The outflow is divided into fast runoff and slow runoff by a distributor according to the fast flow distribution coefficient frac_fast. The fast runoff enters the fast power function reservoir, and the slow runoff enters the slow power function reservoir. The outflow from the two reservoirs is added together to obtain the total runoff.

[0093] Step 34, GR4J structure pre-assembly: interception filter → runoff reservoir → unit line → runoff collection module, French classical hydrological model, containing the complete runoff generation and collection process, with 4 parameters, specifically including runoff reservoir capacity x1, groundwater exchange coefficient x2, runoff collection module capacity x3, unit line lag time lag_time, and its specific connection topology is as follows:

[0094] After pretreatment by the interception filter, the net precipitation and net evapotranspiration enter the runoff reservoir. The infiltration amount of the runoff reservoir is distributed to two unit lines, fast and slow, according to the ratio ψ. The lag time of the fast unit line is τ, and the lag time of the slow unit line is 2τ. The output of the fast unit line enters the runoff module. The outflow of the runoff module is added to the output of the slow unit line to obtain the total runoff.

[0095] Step 35, HYMOD structure pre-assembly: Upper soil zone → water distributor → 3 tandem linear reservoirs (fast path) + 1 linear reservoir (slow path), using a probability distribution to represent spatial heterogeneity, with 5 parameters, specifically including maximum soil water storage capacity Smax, runoff index beta, fast path diversion ratio frac_fast, fast linear reservoir runoff coefficient k_fast, and slow linear reservoir runoff coefficient k_slow. The specific connection topology is as follows:

[0096] Rainfall and potential evapotranspiration enter the upper soil zone. The runoff is divided into fast and slow paths by a distributor according to the fast flow distribution coefficient frac_fast. The fast runoff passes through three linear reservoirs in series, while the slow runoff enters one linear reservoir. The total runoff is obtained by adding the outflows at the end of the fast path and the outflows of the slow path.

[0097] Step 36, M4 structure pre-assembly: unsaturated reservoir → power function reservoir, adding runoff index parameters to the simple HBV structure, 5 parameters, specifically including maximum soil water storage capacity Smax, evapotranspiration coefficient Ce, runoff generation index beta, runoff coefficient k, and power exponent of storage-discharge relationship alpha;

[0098] Step 37: Pre-assembly of the two-layer water tank structure: net precipitation → water tank 1 → water tank 2, with 6 parameters, specifically including the outflow coefficient a11 of the side hole of the first layer water tank, the threshold height h11 of the side hole of the first layer water tank, the outflow coefficient a12 of the bottom hole of the first layer water tank, the outflow coefficient a21 of the side hole of the second layer water tank, the threshold height h21 of the side hole of the second layer water tank, and the outflow coefficient a22 of the bottom hole of the second layer water tank.

[0099] Step 38: Pre-assembly of the three-layer water tank structure: net precipitation → water tank 1 → water tank 2 → water tank 3, with 9 parameters, specifically including the outflow coefficient a11 of the side hole of the first layer water tank, the threshold height h11 of the side hole of the first layer water tank, the outflow coefficient a12 of the bottom hole of the first layer water tank, the outflow coefficient a21 of the side hole of the second layer water tank, the threshold height h21 of the side hole of the second layer water tank, the outflow coefficient a22 of the bottom hole of the second layer water tank, the outflow coefficient a31 of the side hole of the third layer water tank, the threshold height h31 of the side hole of the third layer water tank, and the outflow coefficient a32 of the bottom hole of the third layer water tank.

[0100] Step 39: Pre-assembly of the four-layer water tank structure: net precipitation → water tank 1 → water tank 2 → water tank 3 → water tank 4, with 12 parameters, specifically including the outflow coefficient a11 of the side hole of the first layer water tank, the threshold height h11 of the side hole of the first layer water tank, the outflow coefficient a12 of the bottom hole of the first layer water tank, the outflow coefficient a21 of the side hole of the second layer water tank, the threshold height h21 of the side hole of the second layer water tank, the outflow coefficient a22 of the bottom hole of the second layer water tank, the outflow coefficient a31 of the side hole of the third layer water tank, the threshold height h31 of the side hole of the third layer water tank, the outflow coefficient a32 of the bottom hole of the third layer water tank, the outflow coefficient a41 of the side hole of the fourth layer water tank, the threshold height h41 of the side hole of the fourth layer water tank, and the outflow coefficient a42 of the bottom hole of the fourth layer water tank.

[0101] Step 4, Parameter Space Definition: Based on the characteristics of each model structure, establish a parameter vector and its physical constraint range for each structure. The parameters include common parameters, HBV model structure parameters, GR4J model structure parameters, HYMOD model structure parameters, and TANK model structure parameters.

[0102] Step 41: Set the parameter range for the common parameters, specifically including the maximum soil water storage capacity Smax, evapotranspiration coefficient Ce, runoff index beta, and evapotranspiration smoothing threshold parameter m, where: the parameter range of Smax is [100.0, 400.0], the parameter range of Ce is [0.5, 0.95], the parameter range of beta is [1.0, 4.0], and the parameter range of m is [0.001, 0.1].

[0103] Step 42: Set the parameter ranges for the HBV model structure parameters, specifically including the drainage coefficient k, the drainage coefficient k1 of the first reservoir, the drainage coefficient k2 of the second reservoir, the drainage coefficient k_fast of the fast reservoir, the drainage coefficient k_slow of the slow reservoir, the runoff shape parameter alpha, the power exponent alpha1 of the first reservoir, the power exponent alpha2 of the second reservoir, and the fast path splitting ratio frac_fast, where: the parameter range of k is [0.02, 0.3], the parameter range of k1 is... The parameter ranges are [0.08, 0.4] for k1, [0.005, 0.08] for k2, [0.08, 0.4] for k_fast, [0.005, 0.04] for k_slow, [0.7, 1.5] for alpha, [0.7, 1.5] for alpha1, [0.7, 1.5] for alpha2, and [0.45, 0.75] for frac_fast.

[0104] Step 43: Set the parameter ranges for the GR4J model structure parameters, specifically including the runoff reservoir capacity x1, groundwater exchange coefficient x2, runoff module capacity x3, and unit linear lag time lag_time, where: the parameter range of x1 is [100.0, 800.0], the parameter range of x2 is [-3.0, 3.0], the parameter range of x3 is [20.0, 250.0], and the parameter range of lag_time is [0.8, 4.0].

[0105] Step 44: Set the parameter range for the HYMOD model structure parameters, specifically including the maximum soil water storage capacity Smax, runoff index beta, fast path diversion ratio frac_fast, fast linear reservoir drainage coefficient k_fast, and slow linear reservoir drainage coefficient k_slow, which have the same value range as the previously set Smax, beta, frac_fast, k_fast, and k_slow.

[0106] Step 45: Set the parameter ranges for the water tank model structure parameters, specifically including the outflow coefficients a11, h11, a12, a21, h21, and a22 for the side holes of the first-layer water tank, the threshold height h31, and the bottom hole of the first-layer water tank, respectively; the outflow coefficients a31, a32, and a41 for the side holes of the fourth-layer water tank, and the threshold height h41 and bottom hole of the fourth-layer water tank, respectively. The parameter range for a11 is [0.15, 0.8]. The parameter ranges for a12 are [0.02, 0.4], a21 are [0.08, 0.6], a22 are [0.01, 0.25], a31 are [0.02, 0.4], a32 are [0.005, 0.15], a41 are [0.01, 0.25], a42 are [0.001, 0.08], h11 are [8.0, 40.0], h21 are [8.0, 40.0], h31 are [8.0, 40.0], and h41 are [8.0, 40.0].

[0107] Step 5 involves using genetic programming for cooperative search, specifically including the following steps:

[0108] Step 51, Genetic Programming Population Initialization: Randomly generate 500 individuals within the population, each individual containing a randomly selected model structure type and a parameter vector uniformly sampled in the corresponding parameter space;

[0109] Step 52, Fitness Assessment: For each individual in the population, hydrological simulations are performed using rate-defined precipitation, potential evapotranspiration, and measured runoff data, and the individual fitness is calculated and assessed.

[0110] Step 53, Individual Selection: A tournament selection strategy is adopted, in which four individuals are randomly selected from the population each time for comparison, and the one with the highest fitness is selected as the parent.

[0111] Step 54, Parent Crossover: For the two selected parent individuals, crossover is performed with probability P. c =0.7 Perform structural crossover and parametric crossover to generate two offspring individuals;

[0112] Step 541, Structural Crossover: Swap the model structure types of the two parent generations with a structural crossover probability of 0.5;

[0113] Step 542, Parameter Crossing: For each parameter required by the offspring, randomly select from the two parent parameters at the corresponding position for inheritance. If a parameter does not exist in the parent due to different structure types, it is randomly initialized in the corresponding parameter space.

[0114] Step 55, Individual Mutation: For offspring generated from the crossover, the mutation probability P is... m =0.25 Perform mutation for each parameter;

[0115] Step 551, Random Perturbation Generation: Generate the variation of the parameter by using a Gaussian random perturbation as a function. ,in σ is the coefficient of variation, with a value of 0.5;

[0116] Step 552, Parameter Update and Boundary Exceedance Handling: According to The parameters are updated. If the mutated parameters exceed the boundary, a mirror reflection strategy is used to bring them back to the valid range.

[0117] Step 56, Elite Retention: Sort the individuals in the current population from high to low fitness, retain the top ten elite individuals to directly enter the next generation population to ensure that the optimal solution is not lost, and replace the remaining individuals with the offspring generated by the crossover and mutation operation.

[0118] Step 57, Early Stop Judgment: Record the optimal fitness of each generation. If the optimal fitness of the population does not improve for 15 consecutive generations, the early stop condition is triggered to terminate evolution.

[0119] Step 58, Iterative Update: Repeat steps 52 to 57 until the maximum number of generations of evolution is reached (100 generations) or the early stopping condition is triggered.

[0120] Step 6: Performance Evaluation: Daily runoff simulations are performed using the optimal model structure and parameters determined through search and calibration to evaluate the model's generalization performance. This includes the following steps:

[0121] Step 61: Use the calibration set data divided in Step 1 to perform hydrological simulation on the population generated in Step 5, and set a warm-up period of 30 days to allow the model state to stabilize. The simulation results during the warm-up period are not included in the subsequent index calculation, thereby reducing the impact of the uncertainty of the initial state.

[0122] Step 62: Calculate the Nash efficiency coefficient (NSE) for each individual. The specific calculation formula is shown below:

[0123]

[0124] Among them: Q obs,i and Q sim,i These represent the observed flow rate and the simulated flow rate, respectively; n is the total time step of the calibration set data.

[0125] Step 63: Introduce a complexity penalty term (CP) based on the BIC concept to control model complexity and avoid overfitting. The specific formula is shown below:

[0126]

[0127] Where: λ is the regularization strength, with a value range of [0.03 – 0.10]; n is the total time step of the calibration set data; and k is the number of parameters in the generated model.

[0128] Step 64: Calculate the final model fitness function and use it to evaluate fitness. The specific formula for calculating the fitness function is shown below:

[0129] .

[0130] A 30-day warm-up period was used to calculate the NSE (Neural Sequence of Performance) index of the optimal model on the validation set to evaluate the model's generalization performance. Simultaneously, the indices of the calibration set and validation set were compared to evaluate the model's fit. The optimal model output included: the optimal model structure type, the optimized parameter vector, the performance indices of the calibration set and validation set, and the runoff simulation sequence. The optimal model structure and its specific parameters are as follows: Figure 2 As shown; the observed flow sequence is compared with the flow sequence simulated by the optimal model, and the results are as follows. Figure 3 As shown; the evolution process of the model structure and parameters is as follows Figure 4 As shown. By Figure 2 It can be seen that the optimal model structure is hbv_parallel with 6 parameters. Comparing the performance of the optimal model on the training set and the validation set: the performance NSE on the calibration set is 0.71 and the performance NSE on the validation set is 0.73. Both the NSE on the validation set and the NSE on the calibration set are higher than 0.7, indicating that the optimal model can simulate the rainfall-runoff process in the Red Creek watershed very well.

Claims

1. A method of constructing a rainfall runoff model, characterized by, Includes the following steps: Step 1: Collect and preprocess relevant daily sequence data within the watershed; Step 2: Construct a hydrological model component library; Step 3: Pre-assemble the hydrological model structure; Step 4: Define the parameter space; Step 5: Use genetic programming for cooperative search; Step 6: Use the optimal model structure and parameters obtained through search and calibration to simulate daily runoff, evaluate the model's generalization performance, and output the optimal model structure type, optimized parameter vector, calibration set, validation set performance index, and runoff simulation sequence.

2. The method of claim 1, wherein, In step 1, the daily rainfall sequence of each rain gauge station in the watershed, the daily outlet cross-sectional flow sequence of the watershed, and the daily net radiation sequence, daily soil heat flux sequence, daily average temperature sequence, daily wind speed sequence at 2m height, and daily actual water vapor pressure sequence used for the calculation of potential evapotranspiration in the watershed are collected. The obtained sequences are then preliminarily processed and divided into calibration set and validation set.

3. The method of claim 1, wherein, In step 2, a hydrological model component library is constructed, which includes eight basic hydrological components, including unsaturated reservoirs, power function reservoirs, linear reservoirs, interception filters, runoff-generating reservoirs, confluence modules, unit lines, and tank units. Each component implements independent water balance calculation and status update functions.

4. The method of claim 1, wherein, In step 3, the hydrological model structure is pre-assembled. Based on the hydrological component library, nine pre-assembled model structures are constructed, including HBV simple structure, HBV series structure, HBV parallel structure, GR4J structure, HYMOD structure, M4 structure, two-layer water tank structure, three-layer water tank structure and four-layer water tank structure. Each structure defines a specific component connection topology and data flow direction.

5. The method of claim 1, wherein, In step 4, the parameter space is defined, and a parameter vector and its physical constraint range are established for each structure. The parameters include common parameters, HBV model structural parameters, GR4J model structural parameters, HYMOD model structural parameters, and TANK model structural parameters.

6. The method of claim 1, wherein, Step 5 involves using genetic programming for cooperative search, specifically including the following steps: Step 51, Genetic Programming Population Initialization: Randomly generate individuals within the population. Each individual contains a randomly selected model structure type and a parameter vector uniformly sampled in the corresponding parameter space. Step 52, Fitness Assessment: For each individual in the population, hydrological simulation is performed using precipitation, potential evapotranspiration and measured runoff data within the set, and the individual fitness is calculated and assessed. Step 53, Individual Selection: A tournament selection strategy is adopted, in which a specified number of individuals are randomly selected from the population each time for comparison, and the one with the highest fitness is selected as the parent. Step 54, Parent Crossover: For the two selected parent individuals, perform structural crossover and parameter crossover with predefined crossover probabilities. Structural crossover swaps the model structure types of the two parents with predefined probabilities. For parameter crossover, for each parameter required by the offspring, randomly select and inherit from the parameters of the two parents at the corresponding position. Step 55, Individual Mutation: For offspring individuals generated by crossover, Gaussian mutation is performed on each parameter with a predefined mutation probability. The mutation amount follows a normal distribution with a mean of 0 and a standard deviation of σ×(parameter upper bound - parameter lower bound), where σ is the mutation intensity coefficient. When the parameter exceeds the boundary after mutation, a mirror reflection strategy is used to bring it back to the effective range. Step 56, Elite Retention: Sort the individuals in the current population from high to low fitness, retain a predefined number of elite individuals to directly enter the next generation population, and replace the remaining individuals with the offspring generated by the crossover and mutation operation. Step 57, Early Termination Judgment: Record the optimal fitness of each generation. If the optimal fitness of the population does not improve for a certain number of consecutive generations, then terminate the evolution. Step 58, Iterative Update: Repeat the above steps until the predefined maximum number of generations is reached or the early stopping condition is triggered.

7. The method of claim 1, wherein, Step 6, evaluating the model's generalization performance, specifically includes the following steps: Step 61: Use the rate set to perform hydrological simulation on the generated population, set a simulation warm-up period to allow the model state to stabilize, and the simulation results during the warm-up period will not be included in the subsequent index calculation; Step 62: Calculate the NSE, KGE, RMSE and VE indices for each individual; Step 63: Introduce a complexity penalty term CP based on the BIC concept to control model complexity and avoid overfitting. The specific formula is shown below: Where: λ is the regularization strength; n is the total time step of the calibration set data; k is the number of parameters of the generated model; Step 64: Calculate the final model fitness function and use it to evaluate fitness.

8. A system for implementing the method of constructing a rainfall runoff model as claimed in claim 1, characterized in that, include: Data acquisition module, component construction module, collaborative search module, and model output module; The data acquisition module collects and preprocesses relevant daily sequence data within the watershed; the component construction module builds a hydrological model component library, pre-assembles the hydrological model structure, and defines the parameter space; the collaborative search module uses genetic programming to perform collaborative search. The model output module uses the optimal model structure and parameters obtained through search and calibration to simulate daily runoff, evaluate the model's generalization performance, and output the optimal model structure type, optimized parameter vector, calibration set, validation set performance index, and runoff simulation sequence.

Images