Impact of climate change on monthly streamflow in the Verde River Basin using two hydrological models

This study assessed the impact of climate change on monthly streamflow in the Verde River Basin, located in the Grande River Basin headwater. For this purpose, the SWAT and VIC hydrological models were used to simulate the monthly streamflow under RCP4.5 and RCP8.5 scenarios, obtained by Regional Climate Models Eta-HadGEM2-ES, Eta-CanESM2 and Eta-MIROC5 in the baseline period (1961-2005) and three time-slice (2011-2040, 20412070, and 2071-2099) inputs. At the end of the century, the Eta-HadGEM2-ES showed larger decrease of precipitation in both radiative scenarios, with an annual reduction of 17.4 (RCP4.5) and 32.3% (RCP8.5), while the Eta-CanESM2 indicated major warming, with an annual increase of 4.7 and 10.2°C under RCP4.5 and RCP8.5, respectively. As well as precipitation changes, the Eta-HadGEM2-ES also showed greater impacts on streamflow under RCP4.5 for the first time-slice (2011-2040), with an annual decrease of 58.0% for both hydrological models, and for the RCP8.5 scenario by the end the century (2071-2099), with an annual reduction of 54.0 (VIC model) and 56.8% (SWAT model). Regarding monthly streamflow, the Eta-HadGEM2-ES and Eta-CanESM2 inputs indicated decrease under the RCP8.5 scenario by the end the century, varying from 7.2 to 66.3 % (VIC model) and 37.0 to 64.7% (SWAT model). In general, Eta-MIROC5 presented the opposite in terms of direction in the simulations with both hydrological models at the end of the century. Combined effects of climate models, hydrological model structures and scenarios of climate change should be considered in assessments of uncertainties of climate change impacts.


INTRODUCTION
Climate change and human activity influence hydrological processes, causing critical impacts on both the dynamics of ecosystems and society, and affect urban water supply systems, irrigated agriculture, and hydropower production (Viola et al., 2014;Oliveira et al., 2017;Byun et al., 2019).
The upper Rio Grande Basin, located in southern Minas Gerais state, has great relevance in water availability and hydropower production, contributing also to the Furnas hydropower plant reservoir (Pinto et al., 2019;Oliveira et al., 2017). Thus, the investigation of regional variations caused by climate change are important for water-resource management, agricultural and energy production, water availability, and forecasting of flooding and drought (Byun et al., 2019).
In South America, some studies (Chou et al., 2014;Ribeiro Neto et al., 2016;Oliveira et al., 2017) evaluated projections Eta driven by HadGEM2 and MIROC5, for two scenarios (RCP4.5 and RCP8 .5). These studies indicated an increase in temperature and a reduction in precipitation, mainly in the central and southeastern regions of Brazil, a fact that may indicate a considerable reduction in streamflow and hence in the potential for energy generation, throughout the 21st century. In addition, it was found in these studies that Eta-HadGEM2-ES is more sensitivity in these changes than the Eta-MIROC5. Alvarenga et al. (2018) evaluated the impact of climate change on the hydrologic components in a watershed located in Mantiqueira Range, southern Minas Gerais state. The authors used the Distributed Hydrology Soil Vegetation Model (DHSVM) hydrological model forced by Eta-HadGEM2-ES under RCPs 4.5 and 8.5, finding a reduction in mean seasonal Impact of climate change on monthly streamflow … Rev. Ambient. Água vol. 16 n. 4, e2683 -Taubaté 2021 streamflow for all future periods and both emission radiative scenarios. The authors showed that most critical situations occur during the summer (December to February) from 2011 to 2040under RCP4.5 and 2071to 2099 These studies evaluated the hydrological impacts in Brazilian basins using downscaling projections from different climate models and scenarios; however, they were limited to one hydrological model. Projections based on one hydrological model may be highly uncertain, while application of different models can show the sources of uncertainties. The more comprehensive the study in sampling uncertainties, the higher its quality (Hattermann et al., 2017;Teklesadik et al., 2017;Krysanova et al., 2018).
Hydrological models can show similar performance during calibration; however, the mean discharge response to climate change may vary due to dissimilarities in hydrological model structure and process equations. In addition, there is a cascade of uncertainty in climate change impact assessment that begins with the construction of future climate scenarios and ends in impact assessment (Bastola et al., 2011;Krysanova et al., 2017). Mach et al. (2017) also highlights the attention in the bias correction, because uncertainties are recognized in the downscaling models. Thus, the use of an ensemble of climate models and a variety of scenarios should be considered to represent a range of impacts.
In this context, the objective of this research was to compare climate change responses in the Verde River Basin, southeast Brazil, using two hydrological models: Variable Infiltration Capacity (VIC) that divides the basin into regular grid cells, and Soil and Water Assessment Tool (SWAT) that disaggregate a basin into subbasins and hydrological response units based on land use, topography and soil classes. The projection data are downscaling of three GCMs (CanESM2, HadGEM2 and MIROC5) by the RCM Eta under two future scenarios (RCP4.5 and RCP8.5). The novelty herein is the application of mesoscale (SWAT) and macroscale (VIC) hydrological models for assessment of climate change on a mesoscale basin, and evaluation of the impacts projected by climate models in the Verde River Basin.

Study area
The Verde River Basin (VRB) has a drainage area of approximately 4,100 km 2 , with an elevation ranging from 809 to 2742 m, and is located in the Mantiqueira Range region, the headwater of the Grande River Basin (GRB), southern Minas Gerais state ( Figure 1). The GRB is an important Brazilian river for hydroelectric energy production (Oliveira et al., 2018;Viola et al., 2014). According to the Kӧppen classification, the climate in this region is mostly Cwb (Alvares et al., 2014), the mean annual precipitation is approximately 1500 mm, and the annual mean temperature is 18°C, with a dry winter season (Mello et al., 2012).
The SWAT is a hydrological model of spatial scale that is continuous, physically-based and semi-distributed in space, which divides the basin into sub-basins connected by a stream network. Land cover, soil type and slope information are combined into hydrological response units (HRUs) consisting of lumped land areas (Arnold et al., 1998). The SWAT was developed to estimate the impact of soil management practices on water, sediment and agricultural chemical productions in large complex river basins with different soil types, soil use and management conditions over long periods of time (Neitsch et al., 2005). The SWAT simulations on the hydrology module are based on water balance (Equation 1), calculated using hydrological routines that calculate the components of the water cycle, such as surface and groundwater flows, canopy storage, evapotranspiration, percolation, lateral flows, and transmission losses. Thus, runoff is calculated for each HRU and routed to estimate streamflow in each sub-basin.
Where SWt and SW0 are, respectively, final and initial soil water content, t is the time, R is the precipitation, Qsurf is the surface streamflow, Ea is the evapotranspiration, wseep is the percolation and bypass flow, and Qgw is the groundwater flow.
The surface flow is calculated by either the SCS Curve Number procedure (Mockus, 1964) using daily precipitation data, or by the Green and Ampt infiltration method (Green and Ampt, 1911). Evapotranspiration is estimated by Penman-Montheith (Montheith, 1965), Priestley-Taylor (Priestley and Taylor, 1972) or Hargreaves (Hargreaves et al., 1985) methods, varying the required meteorological inputs for each method. In this study, the SCS curve number and Penman-Montheith methods were adopted. More details about SWAT theoretical equations can be found in Neitsch et al. (2011) and Arnold et al. (1998).
The VIC model is a grid-based macroscale semi-distributed land surface hydrological model, which consists of two modules: (i) rainfall-runoff (Liang et al., 1994(Liang et al., , 1996Liang and Xie, 2001) that simulates the water and energy fluxes that govern the terrestrial hydrological cycle; and (ii) routing model (Lohmann et al., 1996(Lohmann et al., , 1998) that calculates the discharge using linearized de Saint-Venant equations. In the upper soil layer, the variable soil moisture capacity curve (Equation 2) controls evaporation, infiltration, and surface runoff; while the lower layer controls the baseflow generation, using the ARNO formulation (Todini, 1996). The evapotranspiration (vegetation transpiration, canopy evaporation, and evaporation from bare soil) is based on the Penman-Monteith equation.
Where i and im are the point and maximum point soil moisture capacity, respectively; A is the fraction of an area for which the soil moisture capacity is less than or equal to i; and bi is the soil moisture capacity shape parameter.
VIC model was manually calibrated using a spatial resolution of 0.01° with a total of 3728 grid cells and was carried out by adjusting the parameters described: i) infiltration parameter (b_infilt) describes the amount of available infiltration capacity as a function of relative saturated gridcell area. A higher value of b_infilt gives lower infiltration and yields higher surface runoff; ii) thickness of the third soil layer (depth_3), interfering in the transpiration (depending on root depth) and baseflow. Where thick layers of soil have slower runoff responses (predominance of baseflow) and high evapotranspiration, however, the result is high moisture retention and high baseflow in dry periods; iii) maximum baseflow that can occur from the lowest soil layer (Ds); and iv) fraction of the maximum soil moisture (Ws). A higher value of Ws will raise the water content required for rapidly increasing, non-linear baseflow, which will tend to delay runoff peaks. In addition to the calibration of the VIC model, the parameters of the surface-flow propagation model (Route) were also refined (kinematic wave speed (C) and diffusivity (D)) (Alvarenga et al., 2020). While the SWAT model was calibrated using the Sequential Uncertainty Fitting Algorithm (SUFI2) implemented in SWATCUP by Abbaspour et al. (2007), with a total of 57 sub-basins and 1503 HRUs. The calibration was done by adjusting 14 parameters that represent the surface and subsurface flow processes. These parameters are: soil evaporation compensation coefficient (ESCO); initial SCS runoff curve number for moisture condition II (CN2); the baseflow recession constant (ALPHA_BF), Groundwater delay time (GW_DELAY); threshold depth of water in the shallow aquifer required for return flow to occur (GWQMN); maximum canopy storage (CANMX); effective hydraulic conductivity in main channel (CH_K2), Manning's "n" value for the main channel (CH_N2); plant uptake compensation factor (EPCO); groundwater "revap" coefficient (GW_REVAP); threshold depth of water in the shallow aquifer for "revap" or percolation to the deep aquifer to occur (REVAPMN); soil available water capacity (SOL_AWC); saturated hydraulic conductivity (SOL_K) and surface runoff lag coefficient (SURLAG) (Alvarenga et al., 2020).
Information from the ASTER sensor with a resolution of 1 arc-second (~30x30 meters) was used to obtain the Digital Elevation Model (DEM) ( Figure 1); images of the LANDSAT 8 sensor with a resolution of 1 arc-second were used to obtain the land-use map through supervised and object-oriented classification techniques; and soil class mapping was provided by State Foundation for the Environment (FEAM, 2010). The atmospheric variables were obtained from the meteorological stations from National Institute of Meteorology, and observed streamflow from the water level station at Três Corações, obtained from the National Water Agency (Figure 1). The hydrological models are run for the period 1990-2005: first 3-yrs were used for warm up, 1993 to 1999 calibration period, and 2000 to 2005 for validation period.
Four statistical indices were calculated in the performance assessment: (i) coefficient of determination (R 2 ) that is used to describe the variance between simulated and observed streamflow; (ii) Nash-Sutcliffe efficiency coefficient (CNS) that is used to reflect the matching degree between simulated and observed streamflow (Nash and Sutcliffe, 1970); (iii) Nash-Sutcliffe efficiency of the logarithmic values of streamflow (CNS Log) that is used to assess the model performance during recession period; and (iv) Relative error (Pbias) that is employed to measure the mean tendency of the difference between simulated and observed streamflow (Gupta, 1999). According to Moriasi et al. (2007), the R 2 and CNS values above 0.50 are considered acceptable, and Pbias less than |25%| present satisfactory results. CNS Log values from 0.35 to 0.50 and 0.5 to 0.7 indicate medium and good hydrological model performance, respectively.

Climate Models, Downscaling and Bias Correction
To meet the goals of this study, two hydrological models (SWAT and VIC) will be created with results from three climate models and two emission radiative scenarios. The climate simulations used in this work were based on the dynamical downscaling of the CanESM2, HadGEM2-ES and MIROC5 global climate model (GCMs) simulations using the Eta regional climate model Marengo et al., 2012;Mesinger et al., 2012;Pesquero et al., 2010). These 3 models were recommended for the assessment of climate change projections in South America (Dereczynski et al., 2020).
The main characteristics of the GCMs are summarized in Table 1. The chosen emission radiative scenarios were: RCP4.5 that is a medium-to-low scenario, assuming a stabilization of radiative forcing to 4.5 W m -2 by 2100; and RCP8.5 that is the most pessimistic situation, with radiative forcing to 8.5 W m -2 by 2100. The downscaling method provided simulations with spatial resolution of 20 km, covering the following periods: 1961 to 2005 for baseline, and 2011 to 2099 for future scenarios (RCP4.5 and RCP8.5) that were divided into three time-slices (2011-2040, 2041-2070, and 2071-2099). The monthly variables simulated by Eta model used to assess the potential hydrological impacts on the VRB were precipitation, maximum and minimum temperatures, solar radiation, relative humidity, and wind speed.
However, simulations from RCMs, such as Eta, are subjected to systematic biases (Graham et al., 2007;Rodrigues et al., 2020), mainly caused by errors in conceptualization, discretization and spatial average of climate variables within a grid-cell (Teutschbein and Seibert, 2012). Therefore, the RCM variables biases were corrected based on observed meteorological data, wherein the linear scaling method was used for precipitation, maximum and minimum temperatures as proposed by Lenderink et al. (2007), and for wind speed as proposed by Haddeland et al. (2012).

RESULTS AND DISCUSSION
The coupled simulations are referred to as Eta-HadGEM2-ES-ES, Eta-CanESM2, and Eta-MIROC5, when the Eta model is driven by HadGEM2-ES, CanESM2, and MIROC5, respectively. The variations and the changes are based on comparisons between future scenarios (RCP4.5 and RCP8.5) and baseline period for each RCM. In this study, the wet season is from October to March, while the dry season is from April to September.

Models Performance
Hydrological models were calibrated and validated to the monthly observed streamflow in the Verde River Basin by Alvarenga et al. (2020). Figure 2 shows the monthly mean hydrographs for the SWAT and VIC in the calibration and validation periods, indicating better performance for peak streamflow simulated for VIC, while SWAT was able to capture the behavior of the recession period in the VRB. In research in the Grande River Basin headwater, Viola et al. (2014) and Oliveira et al. (2018), using LASH and SWAT models, respectively, also showed underestimation of simulated streamflow during the validation period. These problems tend to get worse with an increasing number of parameter values. Additionally, in the same region, Alvarenga et al. (2016), applying the DHSVM model, presented underestimation during the dry season and overestimation during the wet season.  Table 2 shows statistical indices from the comparison between observed and simulated monthly streamflow by the SWAT and VIC models. According to Moriasi et al. (2007), the statistical indices showed satisfactory performances predicting the monthly streamflow from VRB, with R 2 and CNS values above 0.50, Pbias less than |25%|, and CNS Log values above 0.4. However, the lower CNS Log value showed poor VIC model performance in the recession period.  Figure 3 shows monthly changes of precipitation and maximum and minimum temperatures projected by RCMs (Eta-CanESM2, Eta-HadGEM2-ES and Eta-MIROC5) for three time-slices (2011-2040, 2041-2070, and 2071-2099) and under RCP4.5 and RCP8.5 scenarios, relative to the baseline period (from 1961 to 2005). The projected changes in monthly precipitation differed among the models, both in terms of the magnitude and in the signal of changes (Figure 3). The precipitation simulations show the greatest decrease of the annual rainy for Eta-HadGEM2-ES when compared to Eta-MIROC5 and Eta-CanESM2 under RCP4.5 and RCP8.5 scenarios.

Climate Model Projections
At the beginning of the 21st century , the Eta-HadGEM2-ES showed a decrease of 32.2 and 23.3% in the annual precipitation under RCP4.5 and RCP8.5, respectively; this situation is also observed at the end of the 21st century (2071-2099), indicating a reduction of 17.4 and 32.3% for RCP4.5 and RCP8.5, respectively. In general, the Eta simulations driven by HadGEM2-ES presented greater decreases during the wet season (October to March period) when compared to the baseline period; whereas June and July, for the third time-slice, revealed an increase of 5.5 and 7.3% for RCP4.5, and 16.2 and 70.5% under RCP8.5, respectively.
Overall, the Eta-HadGEM2-ES produces more precipitation than Eta-MIROC5 during the dry season, and less during the wet season. This behavior can explain the minor decrease in the annual precipitation with Eta-MIROC5 when compared to Eta-HadGEM2-ES. These results are similar to those presented by Chou et al. (2014) in the Central-South region of Brazil, by Lyra et al. (2017) in Central America, and by Alvarenga et al. (2018) and Oliveira et al. (2017) in the Grande River Basin headwater.
Considering the Eta-HadGEM2-ES, Oliveira et al. (2017) observed a reduction in mean annual precipitation of 2.2 and 28.1% for RCPs 4.5 and 8.5, respectively; while regarding the Eta-MIROC5 model, there was an attenuation of 7.5 and 1.2% for RCPs 4.5 and 8.5, respectively. Alvarenga et al. (2018) showed a negative trend in the summer season (December to February) in the same scenarios and Eta-HadGEM2-ES simulations toward the end of the century.
These uncertainties and significant differences also may be attributed to the different setup of GCMs and resolution of RCM Eta, generating responses that diverge in both precipitation forecasting and impacts caused by future scenarios. According to Alvarenga et al. (2018) who used the Eta-5km resolution, the high-resolution (5km) of Eta may capture topographic features and extreme values better when compared with low-resolution (20km).
The effect of increased mean maximum and minimum temperature is evident for three RCMs and in all future scenarios. The Eta-CanESM2 presented greater changes in the monthly maximum and minimum temperatures, with an increase of 7.2 and 4.3°C, respectively, under RCP4.5; and an increase of 13.6 and 7.5°C under RCP8.5, respectively. The Eta-MIROC5 and Eta-HadGEM2-ES showed less increase in the values of temperatures, with an increase of 2.2 and 4.8°C in the maximum temperature, and 2.4 and 2.9°C in the minimum temperature under RCP 4.5; and an increase of 4.3 and 9.2°C in the maximum temperature, and 4.6 and 5.2°C in the minimum temperature, under RCP 8.5 (Figure 3). Figure 4 shows mean monthly hydrographs simulated by SWAT and VIC models for RCMs Eta-CanESM2, Eta-HadGEM2-ES and Eta-MIROC5 during baseline period (grey) and three time-slices (2011-2040, 2041-2070 and 2071-2100), under RCP4.5 (blue) and RCP8.5 (red).

Monthly Streamflow Changes
At the end of the century, the SWAT and VIC model impacts presented similar behaviors in the monthly streamflow in the Eta-HadGEM2-ES, with a decrease varying from 31.9 to 47.0% (RCP4.5) and from 37.0 to 64.7% (RCP8.5) in the simulations with SWAT, and ranging from 15.1 to 58.4% (RCP4.5) and from 44.0 to 61.2% (RCP8.5) in the simulations with VIC. Both hydrological models did not show noticeable differences between the monthly streamflow when compared to the baseline period. Regarding 10nual mean streamflow under Eta-HadGEM2-ES, from 2011to 2040, 2041to 2070, and 2071to 2099and 38.3, 41.3 and 54.0% under the RCP8.5 scenario, respectively. The SWAT model also indicated a reduction in these time-slices of 58. 0,35.2 and 40.5% under RCP4.5,respectively;and 42.1,42.1 and 56.8% under RCP8.5,respectively (Figure 5). In agreement with these results, Alvarenga et al. (2018) and Oliveira et al. (2017), also reported that the Eta model driven by HadGEM2-ES indicated greater changes at the beginning of the 21st century under RCP4.5, while RCP8.5 presented greater changes at the end the 21st century.
On the other hand, the SWAT and VIC simulations with Eta-MIROC5 presented a different behavior in relation to Eta-HadGEM2-ES. At the end of the century, considering the SWAT model, the mean annual streamflow decreased 17.5% (RCP 4.5) and increased 10% (RCP 8.5); additionally, the VIC model indicated a decrease of 23.6% under RCP 4.5 and an increase of 15.8% under RCP 8.5. Under the RCP8.5 scenario, the increase took place mainly in the wet season (October to January) in both hydrological models. However, from January through September there was a reduction of streamflow under RCP4.5 in both hydrological models. Regarding Eta-MIROC5, Oliveira et al. (2017), also by the end of the century, showed a reduction and increase of the streamflow simulations under RCP4.5 and RCP8.5 scenarios, respectively, especially in the wet season (October to January).  The largest differences in streamflow simulations were observed with the Eta-CanESM2. The SWAT model showed, in the annual mean streamflow, a decrease of 12.6% (RCP 4.5) and 51.5% (RCP 8.5) at the end of the century. In the same period, the VIC model indicated a reduction of 39.6% under RCP8.5 and an increase of 12.1% under RCP4.5, which occurred between January and August. Figure 5 shows projected changes (ensemble mean) of the mean monthly streamflow for three time-slices (2011-2040, 2041-2070 and 2071-2100) in relationship to the baseline period (from 1961 to 2005).
At the beginning of the century, both hydrological models presented low variability in the streamflow changes in the future scenarios (RCP4.5 and RCP8.5). The changes can be seen clearly at the end of the century for the RCP8.5 scenario, when compared to the RCP4.5 scenario ( Figure 5).
Hydrological model simulations showed little variation between streamflow outputs. The ensembles mean of the SWAT model indicated a reduction in the streamflow ranging from 13.2 to 30.0% and 22.3 to 40.0% under RCP4.5 and RCP8.5, respectively, by the end of the century. The VIC model also showed a decrease in the streamflow varying from 12.0 to 22.2% (RCP4.5) and 7.5 to 31.9% (RCP8.5). The SWAT model indicates noticeable changes, although it showed a smaller amplitude of monthly maximum and minimum values of streamflow, mainly under RCP 4.5 ( Figure 5).
In general, the predicted average streamflow is relatively closer in the hydrological model responses. However, the amplitude of the simulated streamflow is higher in VIC compared to SWAT. These sharp differences in the output hydrological models can be characterized by seasonal variability in rainfall. Mainly in the rainy season, VIC was more sensitive to peak streamflows, while SWAT reduces maximum and minimum streamflows. The main differences in the streamflow simulations during the rainy season are harmful for agriculture and energy generation in the region. The difference in the parameterization and structure of the hydrological models provided different results of impacts on the average monthly streamflow of the Rio Verde Basin, i.e., the model affects the result of the impact analysis.

CONCLUSIONS
SWAT and VIC models were not calibrated using similar methods, but showed similar performance during calibration according to the statistical precision indices. The hydrological models driven by the RCMs Eta-CanESM2, Eta-HadGEM2-ES and Eta-MIROC5 simulated the monthly streamflow in the Verde River Basin under the influence of RCP4.5 and RCP8.5 scenarios.
There is an overall consensus among all three climate models of a projected maximum and minimum temperature increase in the area. The largest increases occur at the end of the century for all models. The monthly changes in maximum and minimum temperature are noticeably larger for the Eta-CanESM2 than for Eta-HadGEM2-ES and Eta-MIROC5 in the third timeslice, under RCP8.5. In contrast, the monthly changes of precipitation are highly dissimilar among climate models.
Regarding monthly streamflow, both hydrological models showed an increase from October to January for the Eta-MIROC5 under RCP8.5, being more significant in November, with an increase of 49.1 and 99.8% for SWAT and VIC models, respectively. The Eta-HadGEM2-ES and Eta-CanESM2 indicated reduction of streamflow in all months and in both hydrological models, except Eta-CanESM2 under RCP4.5 for VIC model simulations that showed an increase ranging from 3.8 to 42.5% between January and August, and a decrease varying from 10.7 to 35.4% between September and December.
This study highlights sources of uncertainty, such as emission/radiative forcing scenarios, downscaling, bias correction, climate modeling, and hydrological modeling. The results suggest that climate models represent the major source of uncertainty in the simulations. However, the water availability reduction results obtained in the hydrological modeling can cause collapse in both hydropower generation and food production in this region and these impact studies of the hydrological response are important for local communities.