Improving Hydrological Simulations with a Dynamic Vegetation Parameter Framework

Abstract

Many hydrological models incorporate vegetation-related parameters to describe hydrological processes more precisely. These parameters should adjust dynamically in response to seasonal changes in vegetation. However, due to limited information or methodological constraints, vegetation-related parameters in hydrological models are often treated as fixed values, which restricts model performance and hinders the accurate representation of hydrological responses to vegetation changes. To address this issue, a vegetation-related dynamic-parameter framework is applied on the Xinanjiang (XAJ) model, which is noted as Eco-XAJ. The dynamic-parameter framework establishes the regression between the Normalized Difference Vegetation Index (NDVI) and the evapotranspiration parameter K. Two routing methods are used in the models, i.e., the unit hydrograph (XAJ-UH and Eco-XAJ-UH) and the Linear Reservoir (XAJ-LR and Eco-XAJ-LR). The original XAJ model and the modified Eco-XAJ model are applied to the Ou River Basin, with detailed comparisons and analyses conducted under various scenarios. The results indicate that the Eco-XAJ model outperforms the original model in long-term discharge simulations, with the NSE increasing from 0.635 of XAJ-UH to 0.647 of Eco-XAJ-UH. The Eco-XAJ model also reduces overestimation and incorrect peak flow simulations during dry seasons, especially in the year 1991. In drought events, the modified model significantly enhances water balance performance. The Eco-XAJ-UH outperforms the XAJ-UH in 9 out of 16 drought events, while the Eco-XAJ-LR outperforms the XAJ-LR in 14 out of 16 drought events. The results demonstrate that the dynamic-parameter model, in regard to vegetation changes, offers more accurate simulations of hydrological processes across different scenarios, and its parameters have reasonable physical interpretations.

1. Introduction

Recent research indicates that since the end of the last century, vegetation ecology has exhibited varying trends across different regions [1]. The types and growth conditions of vegetation can alter the hydrological mechanisms, such as evapotranspiration and interception, and influence the water balance and runoff processes [2,3]. Although vegetation growth exhibits seasonal variation trends, the parameters closely related to vegetation remain at fixed values in most lumped hydrological models. This means that these models cannot account for the response of hydrological processes to vegetation changes. Previous studies have highlighted that using fixed parameters in hydrological models can limit their simulation performance [4,5,6]. This calls for an eco-hydrological model with dynamic vegetation-related parameters that is capable of simulating the impact of vegetation changes on hydrological processes.

To improve the simulation accuracy of hydrological models, many researchers have employed hybrid methods to model watershed hydrological processes dynamically. Espinoza et al. [7] developed a hybrid model by coupling Long Short-Term Memory (LSTM) networks with the Simple Hydrological Model (SHM), showing that the hybrid model outperforms both the SHM and the LSTM individually. These models demonstrate that the improvement in accuracy stems from the dynamic parameters output by the LSTM. Similarly, Feng et al. [8] established a hybrid HBV model, which showed that the dynamic parameters β and γ allowed it to achieve a simulation performance comparable to that of the LSTM. Other researchers have also employed LSTM to incorporate dynamic parameters into hydrological models [9,10]. By leveraging deep learning black-box models such as LSTM, the hybrid models can achieve the dynamic parameterization of hydrological models without consideration of physical processes through loss-function minimization. At the same time, with hydrological models, hybrid models can still simulate reliable hydrological variables, even if these variables have not been trained in the models, which is not achievable with standalone black-box models [11]. While hybrid models significantly enhance the simulation capability of hydrological variables, they also increase model complexity, which is associated with the increasing uncertainty. This uncertainty manifests not only in the model’s internal behavior but also in the parameters, particularly those with unclear physical mechanisms [12].

Another method for achieving model dynamization is to incorporate additional hydrometeorological and ecological data. For instance, the extra vegetation data, such as the Normalized Difference Vegetation Index (NDVI) or Leaf Area Index (LAI), provide a finer description of dynamic vegetation state variables, which indirectly dynamizes the model instead of modifying the fixed parameters directly. Jin et al. [13] modified the SWAT model using the high-resolution remote-sensing data of LAI, which solved the problem of LAI spatial discontinuity in the SWAT model. The Penman–Monteith (PM) formula has long been regarded as the most reliable method for estimating potential evapotranspiration in watersheds [14,15,16,17]. However, the PM formula cannot accurately partition evapotranspiration contributions between bare soil and vegetation. The impact of vegetation changes on watershed evapotranspiration is significant. For instance, during summer, increased root activity may enhance actual evapotranspiration due to enhanced water uptake, while drought can delay the recovery of water use efficiency [8]. To address this, Leuning et al. [18] improved the PM formula by incorporating LAI (PML formula) to separately calculate evapotranspiration contributions from bare soil and vegetation. Bai et al. [19] replaced the empirical evapotranspiration calculation formula with the PML formula in the Hydroinformatic modeling system (HIMS) model and found that the modified model showed a significant improvement in the streamflow and soil moisture simulations as compared to the original model. This is a reliable and interpretable approach, where hydrometeorological and ecological data are input into the model to correct the state variables based on equations derived from physical processes. However, this approach can only correct the state variables of the model and cannot utilize the dynamical parameters in the model. In reality, the vegetation changes in the hydrological processes are also reflected in model parameters, which controls the simulation of canopy interception, evapotranspiration, and rainfall infiltration [20,21,22,23].

Among various vegetation index, the Normalized Difference Vegetation Index (NDVI) effectively reflects the spectral, morphological, and environmental characteristics of vegetation across different types, seasons, and growth conditions [24]. Due to its versatility, the NDVI has been widely adopted by researchers [25,26,27]. In eco-hydrological models, vegetation-related parameters can be inferred by establishing a regression model with the NDVI [28]. In this study, a regression relationship between the NDVI and the model parameters is established to achieve dynamic parameterization, ensuring the interpretability of dynamic parameters. This method not only fully reflects the vegetation changes within the watershed but also fully considers the dynamic changes of model parameters under vegetation changes. By establishing the regression relationships between the NDVI and vegetation-related hydrological model parameters, the method proposed in this study can effectively introduce dynamic parameters to hydrological models. It also provides a tool to explore the impacts of vegetation changes on a broader range of hydrological processes, extending beyond the rainfall–vegetation change-runoff relationship.

2. Study Area and Data

The Ou River is the second largest river in the Zhejiang Province, bounded by 118°45′ E–121°00′ E and 27°28′ N–28°48′ N (Figure 1). It traverses the mountainous region of southern Zhejiang and flows eastward into the East China Sea. The study area is the basin upstream of the Qingtian Hydrological Station, covering a drainage area of 13,735 km². It is characterized as a mountainous river basin in a subtropical monsoon climate zone, featuring significant elevation variations and dense forest cover [29]. The primary vegetation types include evergreen coniferous forests, broad-leaved forests, mixed coniferous and broad-leaved forests, and a substantial proportion of artificial vegetation [30].

The gauged meteorological data for the operation of the hydrological model, including daily rainfall and evapotranspiration, are collected from five meteorological stations located in and around the Ou River basin as shown in Figure 1 (black points). Areal rainfall and evapotranspiration are calculated through the Thiessen polygon method. The observed discharge data for modeling covers the period from 1 January 1990 to 31 December 1995, as provided by the Qingtian Hydrological Station. This study also uses monthly NDVI data from December 1989 to December 1995, which is sourced from the Global Inventory Modeling and Mapping Studies (GIMMS). The GIMMS NDVI dataset has a spatial resolution of 8 km × 8 km and a temporal resolution of 15 days [31]. The NDVI data is aggregated on a monthly scale. In this study, the calibration and validation periods for the hydrological models are 1990–1993 and 1993–1995, respectively.

3. Methods

Figure 2 shows the flowchart of this study. In this study, an eco-hydrological model based on the Xinanjiang (XAJ) model, called Eco-XAJ, is proposed, which incorporates vegetation change into the model using an NDVI regression module. Then, the discharge simulation performance of the Eco-XAJ model and the XAJ model is compared during the full period, the dry seasons, and the drought events in the Ou River basin.

3.1. Hydrological Model

The Xinanjiang (XAJ) model is a lumped hydrological model developed in the 1970s [32]. The XAJ model operates on a daily time step and is based on the principles of the water balance with four modules, including evapotranspiration, runoff generation, runoff separation, and routing. It divides the watershed into several sub-basins or hydrological response units, each of which is represented by a set of parameters that govern its hydrological behavior. The model is especially suitable for catchments in humid and semi-humid areas, where rainfall is often intense and variable. Key features of the XAJ model include its simplicity, its adaptability to different environmental conditions, and its ability to incorporate spatial variations in rainfall and other meteorological variables. It has been widely applied in flood forecasting, water resource management, and environmental studies, making it a critical tool for hydrologists and researchers seeking to understand and manage water systems in regions with diverse hydrological characteristics [33,34].

The unit hydrograph method is used as the routing modules in the XAJ model [32]. This method applies the unit hydrograph for surface runoff, while interflow and groundwater flow are routed through linear reservoirs. After routing to the outlet section, the results of the linear reservoir routing are combined with the unit hydrograph runoff, and the final simulated discharge at the catchment outlet is obtained. Another method is the Linear Reservoir method, where the surface runoff, interflow, and groundwater flow are all routed through linear reservoirs [33]. In order to determine the suitable runoff routing method for the study area, this study applies both the unit hydrograph method and the Linear Reservoir method to route the runoff generated by the XAJ model, referred to as XAJ-UH and XAJ-LR, respectively.

3.2. Dynamic Parameterization Considering Vegetation Changes

In the XAJ model, the conceptual parameter K directly affects the model’s evapotranspiration, thereby influencing subsequent calculations of soil moisture and runoff generation. As shown in Equation (1), K is a correction factor from potential evapotranspiration to actual evapotranspiration:

$$\mathrm{E}=\mathrm{K}\times {\mathrm{E}}_{\mathrm{p}}$$

(1)

where $\mathrm{E}$ and ${\mathrm{E}}_{\mathrm{p}}$ represent the actual evapotranspiration and potential evapotranspiration, respectively. Considering the impact of vegetation on evapotranspiration [35], the modified XAJ model assumes that the value of K is dynamic and related to the NDVI, in contrast to the original model where K is assumed to be a constant. Here, a regression relationship between K and the NDVI is established and then integrated into the XAJ model, which is referred to as the Eco-XAJ model.

3.3. Calibration and Evaluation of the Model

The parameters of the XAJ model are calibrated with the Non-Dominated Sorted Genetic Algorithm II (NSGA-II) [36,37]. A more detailed description of the model’s calibration and sensitivity analysis has been provided by Liu et al. [33], and the basic information of parameters calibrated in our study is listed in

Table A1. Parameters to be calibrated and sensitive parameters (marked with bold) of the XAJ model.

Module	Parameter	Physical Meaning	Range
Evapotranspiration	K	Ratio of potential evapotranspiration to pan evaporation	0.7–1.1
	C	Evapotranspiration coefficient of the deeper soil layer	0.01–0.05
Runoff generation	WUM	Soil water tension capacity of the upper layer	0–100
	WLM	Soil water tension capacity of the lower layer	0–100
	WDM	Soil water tension capacity of the deep layer	0–50
	B	Exponential of the distribution to water tension capacity	0.1–0.8
	IMP	Percentage of impervious and saturated areas in the catchment	0–0.05
Runoff separation	SM	Free water capacity of the surface soil layer	0–100
	EX	Exponent of the free water capacity curve influencing the development of the saturated area	0–2.5
Routing	KI	Outflow coefficients of soil’s free water storage to interflow	0–0.7
	KG	Outflow coefficients of soil’s free water storage to groundwater	0–0.7
	CI	Recession constants of the interflow	0–1
	CG	Recession constants of the groundwater storage	0–1
Routing (UH)	UH	Unit hydrograph	0–50
	L	Lag time	Defined
Routing (LR)	CS	Recession constants of the surface flow	0–1

.

Three objective functions are used in the model calibration, i.e., the Nash–Sutcliffe Efficiency coefficient ($NSE$), the Correlation Coefficient ($r$), and the Percent Bias ($pBIAS$).

$$NSE=1-{\displaystyle \frac{\sum _{i=1}^N{\left(y_i-{\widehat{y}}_i\right)}^2}{\sum _{i=1}^N{[y_i-\overline{y}]}^2}}$$

(2)

$$r={\displaystyle \frac{\sum _{i=1}^N\left[{\widehat{y}}_i-\stackrel{̿}{y}\right]\left[y_i-\overline{y}\right]}{\sqrt{\sum _{i=1}^N{\left[{\widehat{y}}_i-\stackrel{̿}{y}\right]}^2\sum _{i=1}^N{\left[y_i-\overline{y}\right]}^2}}}$$

(3)

$$pBIAS={\displaystyle \frac{\sum _{i=1}^N{{\widehat{y}}_i-y}_i}{\sum _{i=1}^N{y}_i}}\times 100\%$$

(4)

where $y_i$ is the observed value, ${\widehat{y}}_i$ is the simulated value, $\overline{y}$ is the mean value of observation, and $\stackrel{̿}{y}$ is the mean value of simulation. According to the definition, when $NSE$ and $r$ approach 1 and $pBIAS$ approaches 0, it indicates that the model simulation is more accurate. In this study, $NSE$ and $r$ are used to evaluate simulation performance. $NSE$ also serves as the basis for filtering parameter K. $r$ is employed to explain the linear correlation between the monthly optimal parameter K and the NDVI, as well as between the dynamic parameter K and the monthly optimal parameter K. $pBIAS$ is used to assess the model’s performance in simulating the water balance.

For the Eco-XAJ model, a two-step calibration method is used. The first step involves calibrating the static model parameters, which is consistent with the XAJ model. The second step is to calibrate the dynamic parameter K. Figure 3 shows the procedure for building the Eco-XAJ model. The specific steps are as follows:

First, calculate the optimal value for each month. The optimal value of K for each month is determined through a traversal with the range set from 0.01 to 7. The remaining model parameters are set as the parameters calibrated in the original model. The parameter K corresponding to the best simulation performance (i.e., the maximum NSE) is selected for each month.

Second, select the effective optimal K. When the XAJ model performs poorly in the simulation due to the limitation of the model structure or because of forcing the data quality, there may be months in which no reasonable simulation results can be obtained, regardless of the changes in parameter K. The optimal K values for these months are treated as anomalies. To avoid the influence of these anomalies, they are excluded based on the following three criteria: (1) the optimal K value for a month falls at the boundary of the specified range; (2) the observed discharge during a month does not align with rainfall patterns; and (3) there is a poor simulation performance during a month, characterized by a substantial deviation between the simulated and observed discharges.

Third, build the dynamic K module. The dynamic K module is built based on the regression relationship between the monthly optimal K parameters and the monthly NDVI sequence. This study also considers the potential lagged effects of the vegetation changes on watershed hydrological processes. Therefore, there are two fitting approaches: (1) fitting the optimal K parameter for the current month with the NDVI of the current month (Lag-0 model) and (2) fitting the optimal K parameter for the current month with the NDVI of the previous month (Lag-1 model). The regression equations derived from the fitting process can calculate the dynamic parameter K, taking the vegetation changes into consideration. Subsequently, replacing the fixed parameter K in the XAJ model with this dynamic parameter K completes the modeling process of Eco-XAJ.

Finally, construct the Eco-XAJ model. The fixed K in the XAJ model is replaced with the NDVI to solve K regression equation in the XAJ model.

4. Results and Discussion

4.1. The Relationship Between K and the NDVI

Figure 4 shows the relationship between the NDVI and the optimal values of K for the same month, as well as the relationship between the current month’s K and the previous month’s NDVI. The calibration and validation results for the XAJ models can be found in

Table A2. Evaluations of the simulation with the parameter sets calibrated in the XAJ model.

Period	Metric	Model
		XAJ-UH	XAJ-LR
Calibration	NSE	0.643	0.569
	r	0.809	0.767
	$pBIAS$ %	−2.2	11.6
Validation	NSE	0.623	0.616
	r	0.815	0.798
	$pBIAS$ %	−21.6	−8.9

. According to Figure 4, for both XAJ-UH and XAJ-LR, K and the NDVI exhibit an exponential relationship. However, the XAJ-LR model shows a more pronounced exponential relationship between K and the NDVI. Additionally, compared to XAJ-UH, XAJ-LR generally shows a stronger correlation between K and the NDVI within the same month. Furthermore, the distribution range of K in the XAJ-LR is broader as compared to the XAJ-UH, with fewer values concentrated between 0.01 and 1. Compared with the XAJ-UH model, the K values in the XAJ-LR model are more concentrated around the upper boundary of 7 rather than the lower boundary of 0.01.

For the XAJ-UH, the regression equation for estimating K in the same month using the same month’s NDVI is shown in Equation (5). The regression equation for estimating the current month’s K using the previous month’s NDVI is shown in Equation (6). For the XAJ-LR, the regression equation for estimating K in the same month using the same month’s NDVI is shown in Equation (7). The regression equation for estimating the current month’s K using the previous month’s NDVI is shown in Equation (8).

$${\mathrm{K}}_{\mathrm{t}}=0.1133{\mathrm{e}}^{2.093{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}_{\mathrm{t}}}+1.134\times {10}^{-15}{\mathrm{e}}^{46.34{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}_{\mathrm{t}}}$$

(5)

$${\mathrm{K}}_{\mathrm{t}}=2.162\times {10}^{-3}{\mathrm{e}}^{8.544\times {\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}_{\mathrm{t}-1}}$$

(6)

$${\mathrm{K}}_{\mathrm{t}}=0.455{\mathrm{e}}^{0.526{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}_{\mathrm{t}}}+8.996\times {10}^{-15}{\mathrm{e}}^{45.17{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}_{\mathrm{t}}}$$

(7)

$${\mathrm{K}}_{\mathrm{t}}=9.17\times {10}^{-3}{\mathrm{e}}^{7.305{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}_{\mathrm{t}-1}}$$

(8)

As indicated by the $r$ values in Figure 4a,c, it is evident that, either using XAJ-UH or XAJ-LR, the fitted curve based on the K and the NDVI of the same month demonstrates a stronger correlation. It makes sense that K increases with the rise in the NDVI. When the NDVI is low, the vegetation in the watershed is not abundant, so the amount of water absorbed by the roots and transpired by the vegetation is also low. Conversely, high values of the NDVI indicate dense vegetation in the watershed, and a large amount of water transpired by the vegetation. The result is consistent with several previous research findings, which indicate that evapotranspiration increases with vegetation greening [38].

Additionally, Figure 4a,c exhibits a shared trend: the fitted curves show an inflection point around an NDVI of approximately 0.68, where the K values start to rise sharply. In contrast, the fitted curves in Figure 4b,d do not display such a pronounced inflection point.

4.2. Simulation Performance of Full Period

The four regression equations are used to build the Eco-XAJ model, where Equations (5) and (6) are used in XAJ-UH, noted as Eco-XAJ-UH with the Lag-0 module and Eco-XAJ-UH with the Lag-1 module, and Equations (7) and (8) are used in XAJ-LR, noted as Eco-XAJ-LR with the Lag-0 module and Eco-XAJ-LR with the Lag-1 module. The simulation results of the Eco-XAJ models and the original XAJ models from January 1990 to December 1995 are compared in

Table 1. Comparison between the Eco-XAJ model and the XAJ model.

Model	NSE	r
XAJ-UH	0.635	0.799
Eco-XAH-UH with Lag-0 module	0.647	0.815
Eco-XAH-UH with Lag-1 module	0.639	0.815
XAJ-LR	0.589	0.768
Eco-XAH-LR with Lag-0 module	0.599	0.777
Eco-XAH-LR with Lag-1 module	0.565	0.767

Note: The best-performing model is highlighted in bold.

. In Table 1, the best-performing Eco-XAJ model is highlighted in bold. According to Table 1, the Eco-XAJ-UH models with both the Lag-0 and Lag-1 modules outperform the XAJ-UH model, while the Eco-XAJ-LR model with the Lag-0 module performs better than the XAJ-LR model. The Eco-XAJ-UH with the Lag-0 module shows the best performance of all six models. Three of the four Eco-XAJ models perform better than the original models. Because the Eco-XAJ models with the Lag-0 module outperform those with the Lag-1 module, the subsequent analysis focuses on the Eco-XAJ models with the Lag-0 module.

Figure 5 and Figure 6 present the hydrographs simulated by Eco-XAJ-UH and Eco-XAJ-LR with the Lag-0 module. Three representative periods are selected, i.e., early May 1990 to early September 1990 (Duration I), early July 1991 to late December 1991 (Duration II), and early July 1995 to late December 1995 (Duration III). The periods of Duration I to III correspond to the wet season, the dry season and the recession periods, respectively. In Duration I, both Eco-XAJ-UH and Eco-XAJ-LR underestimate the flood peaks. The Eco-XAJ-LR model underestimates floods more significantly than the Eco-XAJ-UH model. As shown in Figure 6b, the three flood peaks exceeding 2000 m³/s, simulated by the Eco-XAJ-LR model before 1 July 1990, are all lower than those simulated by the XAJ-LR model. However, when the observed discharge is below 2000 m³/s, the simulations of the Eco-XAJ models are generally better than those of the original XAJ models. Specifically, the XAJ-LR model incorrectly simulates the flood peak on 3 July 1990, while the Eco-XAJ-LR simulates the discharge closer to the observation.

According to Figure 5c and Figure 6c, the XAJ models overestimate the flood peaks in Duration II. The Eco-XAJ models address this issue, mainly due to their ability to simulate more accurate evapotranspiration. According to Figure 5d and Figure 6d, there are several short rainfall events within the recession process in Duration III. It can be observed that the discharge does not show a significant flood peak despite the rainfall reaching its peak. The XAJ models present incorrect flow peaks for all these rainfall events, while the Eco-XAJ models almost perfectly reproduce the variations of discharge during this period, with only one incorrect flow peak. For the rainfall events less than 100 mm, the XAJ models and the Eco-XAJ models incorrectly simulate the flood peaks, but the flood peaks simulated by Eco-XAJ are lower than those of the original model.

To better visualize the differences between the models, the simulation results of the Eco-XAJ and XAJ models are presented in Figure 7. In Figure 7a, the Eco-XAJ model’s scatter points are more concentrated around the 1:1 line, indicating that the simulated discharge of the Eco-XAJ-UH model is closer to the observed discharge, as compared to the original model, when the observed discharge is less than 2000 m³/s. In Figure 7a, it can be observed that the outliers primarily occur when the discharge exceeds 2000 m³/s, where the XAJ model underestimates the flood peak. This is consistent with the results shown in Figure 5a. In Figure 7b, it can be observed that the maximum deviation of discharge between the Eco-XAJ-UH model and the original model is approximately −1000 m³/s, when the discharge is less than 2000 m³/s. As the simulated discharge of the original model exceeds 2000 m³/s, the tendency for the Eco-XAJ-UH model to simulate lower discharges, as compared to the original model, gradually diminishes. When the simulated discharge of the original model surpasses 4000 m³/s, the simulated discharge of the Eco-XAJ-UH model becomes nearly identical to that of the original model. This indicates that the Eco-XAJ-UH model almost does not reduce the original model’s flood-peak simulation performance. The comparison between the Eco-XAJ-LR model and the XAJ-LR model is shown in Figure 7c,d. The conclusions drawn with Eco-XAJ-LR model are similar to those drawn with the Eco-XAJ-UH model. According to Figure 7b,d, when the discharge is below 3000 m³/s, the difference between the Eco-XAJ-LR and the XAJ-LR models is more significant than that between the Eco-XAJ-UH and XAJ-UH models, with the maximum deviation being around −2250 m³/s. This deviation mainly occurs during the wet seasons, as shown in Figure 6b. In Duration I, the XAJ-LR model overestimates the discharge, while the simulation of Eco-XAJ-LR model closely matches the observation. This indicates that the improvement in the evapotranspiration process in the Eco-XAJ model effectively enhances the model’s simulation performance.

For a more detailed analysis, the error plots of the Eco-XAJ model corresponding to the observed discharges are presented in Figure 8. Figure 8a,c illustrates the fact that when the observed discharge is below approximately 1000 m³/s, the absolute error of the Eco-XAJ model is generally lower than that of the XAJ model. As shown in Figure 8b,d, the error analysis of the discharge less than 1000 m³/s is conducted, where a ∆Error greater than 0 indicates that the Eco-XAJ model outperforms the XAJ model. The red solid line shows the cumulative ∆Error. If the slope of the cumulative ∆Error is greater than 0, it indicates that the advantage of the Eco-XAJ model is increasing. According to Figure 8b,d, it is evident that the difference between the Eco-XAJ-LR and XAJ-LR models is more pronounced as compared to that between Eco-XAJ-UH and XAJ-UH. The cumulative ∆Error in Figure 8d is approximately twice that shown in Figure 8b. This indicates that the improvement in the simulation accuracy of Eco-XAJ-LR is greater than that of Eco-XAJ-UH. Furthermore, both subplots reveal that the advantage of the Eco-XAJ models increases with discharge ranging from 0 m³/s to 600 m³/s. However, in the range of the discharge from 600 m³/s to 1000 m³/s, the slope is approximately 0, suggesting that the advantage of the Eco-XAJ models weakens.

4.3. Simulation Performance During Dry Seasons

The hydrographs during the dry seasons (from October to February) are shown in Figure 9 and Figure 10, in which the advantages and limitations of the Eco-XAJ models are demonstrated. The Eco-XAJ model also exhibits notable improvements, especially in correcting the XAJ model’s overestimation of peak discharges. For example, the XAJ models significantly overestimate the peak on 5 October 1991, 10 November 1993, and 12 October 1994, but the Eco-XAJ models successfully mitigates this overestimation, and the Eco-XAJ-LR model shows the most significant improvement and provides the best match for the observed discharge. When the XAJ models fail to simulate the discharge, the improvement of the Eco-XAJ models is also limited. In the dry seasons of 1990 and 1992, both the original XAJ models and the Eco-XAJ models are not able to simulate the peak discharge accurately. In the dry season, the Eco-XAJ models’ runoff simulation performances generally improve. However, in a few cases, due to the poor performance of the original XAJ models, Eco-XAJ does not show significant improvements in the discharge simulation.

To quantitatively analyze the simulation performance of the models during the dry seasons, the $NSE$ and $pBIAS$ values for all models are listed in

Table 2. $NSE$ and $pBIAS$ during the dry seasons.

Metric	Model	Period
		1990	1991	1992	1993	1994
$NSE$	Eco-XAJ-UH	−0.387	0.616	−1.027	−0.054	0.726
	Eco-XAJ-LR	0.225	0.511	−0.222	−0.072	0.682
	XAJ-UH	0.213	−0.020	−1.033	−0.210	0.703
	XAJ-LR	0.343	0.058	−0.669	−0.310	0.654
$pBIAS$ %	Eco-XAJ-UH	−80.3	−14.2	−63.0	−43.5	−28.3
	Eco-XAJ-LR	−3.5	35.1	42.4	32.4	6.9
	XAJ-UH	−61.4	20.2	−63.9	−21.9	−19.8
	XAJ-LR	10.8	61.6	55.8	43.2	12.6

Note: The best-performing model is highlighted in bold.

. It is observed that all four models have negative $NSE$ values and large $pBIAS$ values during the dry seasons in 1992, suggesting a poor performance of the XAJ and Eco-XAJ models, while all models perform well in 1994. The Eco-XAJ-UH model demonstrates the best $NSE$ among the four models during the dry seasons in 1991, 1993, and 1994. Notably, in 1991, Eco-XAJ-UH also achieves the best $pBIAS$. As analyzed earlier, the Eco-XAJ-UH model exhibits superior peak simulation performance as compared to XAJ-UH, avoiding erroneous peaks that do not align with the observed discharge (around 1 November 1991), which contributes to its higher NSE. Furthermore, during the period between 5 October 1991 and 6 January 1992, the Eco-XAJ-UH model simulates a discharge that is closely aligned to the observation, leading to its top $pBIAS$ performance. In 1994, the NSE for Eco-XAJ-UH improves slightly by 0.023 as compared to XAJ-UH, reflecting its success in addressing peak overestimation. The Eco-XAJ-LR model achieves the best $pBIAS$ among the four models during the dry seasons in 1990, 1992, and 1994. However, as previously noted, the simulated hydrographs of all models in 1990 and 1992 are not reasonable. In 1994, similar to the Eco-XAJ-UH model, the Eco-XAJ-LR model shows a slight NSE improvement (0.028) as compared to XAJ-LR, which is also due to its success in addressing peak overestimation, resulting in its optimal $pBIAS$. In 1991, the Eco-XAJ models significantly surpass the XAJ models in both NSE and $pBIAS$ due to the successful correction of the overestimation of peak discharge.

4.4. Simulation Performance of Drought Events

Considering that water balance is more critical during drought events, $pBIAS$ is used as the evaluation metric for discharge in drought events. In this study, drought events are defined as events with a peak discharge of 400 ± 80 m³/s during a 20-day period.

Table 3. $pBIAS$ Comparison of the Eco-XAJ models and the original model during low-flow rainfall events.

Durations	$\mathit{p}\mathit{B}\mathit{I}\mathit{A}\mathit{S}$%				NDVI
	Eco-XAJ-UH	Eco-XAJ-LR	XAJ-UH	XAJ-LR
1990-11-17–1990-12-07	−87.0	−8.7	−53.2	11.5	0.7202, 0.6895
1991-06-25–1991-07-15	94.4	143.3	132.1	173.7	0.6700, 0.7388
1991-07-15–1991-08-04	−72.5	19.4	72.5	90.3	0.7388, 0.6931
1991-08-24–1991-09-13	−44.4	25.2	18.6	56.5	0.6931, 0.7134
1991-10-03–1991-10-23	17.0	38.4	172.5	176.6	0.7483
1992-01-11–1992-01-31	−58.5	75.5	−56.2	90.3	0.6334
1993-02-14–1993-03-06	−25.3	−0.9	−44.7	3.9	0.4822, 0.5836
1993-08-13–1993-09-02	−91.6	−15.9	−67.1	1.5	0.7179, 0.6819
1993-11-21–1993-12-11	−45.6	31.0	−16.3	37.5	0.6328, 0.6613
1994-07-19–1994-08-08	−0.2	47.0	47.3	78.1	0.6877, 0.6755
1994-12-26–1995-01-15	−43.7	−5.1	−45.9	−4.2	0.5556, 0.6637
1995-02-04–1995-02-24	−20.9	−0.8	−7.2	2.4	0.6505
1995-07-14–1995-08-03	−29.0	−10.9	91.3	82.1	0.7303, 0.6899
1995-08-03–1995-08-23	24.6	47.3	126.1	141.3	0.6899
1995-08-23–1995-09-12	23.5	51.3	46.7	94.9	0.6899, 0.7171
1995-09-12–1995-10-02	−70.2	34.9	−57.8	65.7	0.7171, 0.7425

Note: The best-performing model is highlighted in bold.

lists the $pBIAS$ for Eco-XAJ and XAJ in 16 drought events, along with the NDVI values. In Table 3, the better $pBIAS$ values for specific events are highlighted in bold. The Eco-XAJ-UH outperforms the XAJ-UH in 9 out of 16 events, while the Eco-XAJ-LR outperforms the XAJ-LR in 14 out of 16 events. This indicates that dynamic parameterization can improve the model’s performance in simulating discharge in drought events. Moreover, it is also observed that Eco-XAJ-LR has more advantages than Eco-XAJ-UH in simulating discharge in drought events. Table 3 also shows that $∆pBIAS$ ranges from 0 to −33.8% in the seven drought events, whereas the Eco-XAJ-UH model does not show improvement. This means that although Eco-XAJ-UH does not always outperform XAJ-UH in water-balance simulation, the performance gap is relatively small.

Figure 11 illustrates the relationship between the NDVI and the model’s performance during drought events. As shown in Figure 11a, the $pBIAS$ values of the XAJ-UH model are negative when the NDVI is below 0.68 but become mainly positive when the NDVI exceeds 0.68. In contrast, the $pBIAS$ values of the Eco-XAJ-UH model are mostly negative, with positive $pBIAS$ values only occurring when the NDVI is greater than 0.68. Figure 11b shows that the $pBIAS$ values for both the XAJ-LR and Eco-XAJ-LR models are generally positive. When the NDVI is below 0.68, the $pBIAS$ values of the two models are similar. However, when the NDVI surpasses 0.68, a significant difference emerges between the models, with the $pBIAS$ of the Eco-XAJ-LR model moving closer to zero. Figure 11a,b reveals that the XAJ-UH and XAJ-LR models produce relatively large $pBIAS$ values when the NDVI exceeds 0.68. This suggests that the XAJ models face challenges in simulating the water balance accurately during drought events with high NDVIs.

Figure 11c,d, indicates that the Eco-XAJ model outperforms the XAJ model in a simulation in which the value on the vertical axis is greater than 0. Conversely, if the value is less than or equal to 0, it suggests that there is no improvement. As shown in Figure 11c, in cases where the model’s performance improves, the magnitude of the improvement is greater than the magnitude of the decline observed in cases where the performance does not improve. It suggests that the Eco-XAJ-UH model generally performs better than the XAJ-UH model. Furthermore, Figure 11c illustrates that the improvement in the water-balance simulation performance of Eco-XAJ-UH, as compared to XAJ-UH, is relative to the NDVI, with more noticeable improvement when the NDVI exceeds 0.68. In Figure 11d, the Eco-XAJ-LR model shows an improvement in water balance simulation as compared to the XAJ-LR model, with only two events having negative values of $∆pBIAS$ near 0. Moreover, similar to the Eco-XAJ-UH model, when the NDVI exceeds 0.68, the Eco-XAJ-LR model demonstrates a noticeable enhancement in water balance simulation as compared to the XAJ-LR model.

In Figure 11, it is evident that the performance improvement from the XAJ model to the Eco-XAJ model is related to the NDVI. When the NDVI is below 0.68, the XAJ-UH model underestimates the discharge, which implies that the XAJ-UH model overestimates evapotranspiration. Correspondingly, the XAJ-UH model underestimates evapotranspiration when the NDVI is above 0.68. Meanwhile, the XAJ-LR model generally overestimates discharge and underestimates evapotranspiration. The deviation becomes more pronounced when the NDVI is above the threshold of 0.68. For the two Eco-XAJ models, the pronounced performance improvement, as compared to the XAJ models, occurs when the NDVI exceeds 0.68. This threshold is consistent with the inflection point observed in Figure 4a,c. This indicates that the Eco-XAJ models, by dynamically adjusting parameter K, simulate more reasonable evapotranspiration, leading to more accurate simulated discharge.

5. Conclusions

This study developed the eco-hydrological model Eco-XAJ to account for the vegetation changes in the Ou River basin. In the Eco-XAJ model, a regression relationship between the parameter K in the XAJ model and the NDVI was established. This method dynamically adjusts the model parameters to account for vegetation changes. The main conclusion can be drawn as follows:

(1) The simulation accuracy of the Eco-XAJ model is generally better than that of the original model. The differences between the two models are mainly observed in low-discharge simulations, with relatively minor differences in high-discharge simulations. The Eco-XAJ model improves the accuracy in the simulation of discharge below 1000 m³/s.

(2) In the dry season, the advantage of the Eco-XAJ model lies in its ability to overcome the peak overestimation and erroneous peak simulations present in the XAJ model. The advantage in overcoming peak overestimation is most evident in the simulation of unimodal peaks, while the advantage in correcting erroneous peak simulations is reflected in the recession phase after flood events. When the XAJ model produces unreasonable simulations throughout the dry season, the improvements brought by the Eco-XAJ model are limited and not significant. In contrast, when the results of the XAJ model are relatively reasonable during the dry season, the overall performance of the Eco-XAJ model shows slight improvement.

(3) During drought events, the $pBIAS$ of the Eco-XAJ model is typically lower than that of the original model. This advantage is more pronounced in the Eco-XAJ-LR model. Considering vegetation changes, the model’s performance has improved from the perspective of water balance, especially for discharge simulation of drought events with a high NDVI.

It is worth noting that in regions affected by seasonal fires or intensive agriculture, vegetation may undergo rapid changes in both spatial extent and composition. The monthly NDVI may not respond quickly enough to these kinds of land cover changes. As a result, this model may not be able to capture this change. In further studies, to settle the problem, daily vegetation data such as Enhanced Vegetation Index (EVI) or LAI can be introduced into the model to capture more detailed vegetation changes. To incorporate the frequency and intensity of agricultural activities or fires into the eco-hydrological model can also help get rid of this limitation.