# Hydrological Alteration Index as an Indicator of the Calibration Complexity of Water Quantity and Quality Modeling in the Context of Global Change

Roxelane Cakir, Mélanie Raimonet, Sabine Sauvage, Javier Paredes-Arquiola, Youen Grusson, Laure Roset, Maite Meaurio, Enrique Navarro, Miguel Sevilla-Callejo, Juan-Luis Lechuga-Crespo, et al.

<figure><img src="image_1.png" /></figure>

## ► To cite this version:

Roxelane Cakir, Mélanie Raimonet, Sabine Sauvage, Javier Paredes-Arquiola, Youen Grusson, et al. Hydrological Alteration Index as an Indicator of the Calibration Complexity of Water Quantity and Quality Modeling in the Context of Global Change. Water, 2020, 12 (1), pp.115. 10.3390/w12010115. hal-03014343

HAL Id: hal-03014343

https://hal.science/hal-03014343

Submitted on 20 Nov 2020

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Article

# Hydrological Alteration Index as an Indicator of the Calibration Complexity of Water Quantity and Quality Modeling in the Context of Global Change

Roxelane Cakir <fnref n="1" />,* D, Mélanie Raimonet <fnref n="1" />, Sabine Sauvage <fnref n="1" />,* D, Javier Paredes-Arquiola <fnref n="2" />, Youen Grusson <fnref n="3" />, Laure Roset <fnref n="4" />, Maite Meaurio <fnref n="5" />, Enrique Navarro <fnref n="6" />, Miguel Sevilla-Callejo <fnref n="6" />, Juan Luis Lechuga-Crespo <fnref n="5" /> D, Juan Jesús Gomiz Pascual <fnref n="7" /> D, Jose María Bodoque <fnref n="8" /> D and José Miguel Sánchez-Pérez <fnref n="1" /> D

<figure><img src="image_3.png" /></figure>

Received: 15 October 2019; Accepted: 26 December 2019; Published: 30 December 2019

**Abstract:** Modeling is a useful way to understand human and climate change impacts on the water resources of agricultural watersheds. Calibration and validation methodologies are crucial in forecasting assessments. This study explores the best calibration methodology depending on the level of hydrological alteration due to human-derived stressors. The Soil and Water Assessment Tool (SWAT) model is used to evaluate hydrology in South-West Europe in a context of intensive agriculture and water scarcity. The Index of Hydrological Alteration (IHA) is calculated using discharge observation data. A comparison of two SWAT calibration methodologies are done; a conventional calibration (CC) based on recorded in-stream water quality and quantity and an additional calibration (AC) adding crop managements practices. Even if the water quality and quantity trends are similar between CC and AC, water balance, irrigation and crop yields are different. In the context of rainfall decrease, water yield decreases in both CC and AC, while crop productions present opposite trends (+33% in CC and -31% in AC). Hydrological performance between CC and AC is correlated to IHA: When the level of IHA is under 80%, AC methodology is necessary. The combination of both calibrations appears essential to better constrain the model and to forecast the impact of climate change or anthropogenic influences on water resources.

**Keywords:** SWAT; calibration; large-scale; water resource; quality; global change

<footnote n="1">¹ ECOLAB, Université de Toulouse, CNRS, INPT, UPS, 31326 Toulouse, France; melanie.raimonet@gmail.com (M.R.); jose-miguel.sanchez-perez@univ-tlse3.fr (J.M.S.-P.)</footnote>

<footnote n="2">² Research Institute of Water Environmental Engineering (IIAMA), Universitat Politècnica de Valencia, 46022 Valencia, Italy; jparedea@hma.upv.es</footnote>

<footnote n="3">³ Department of Soil and Environment, Faculty of Natural Resources and Agricultural Sciences, Swedish University of Agricultural Sciences (SLU), 750 07 Uppsala, Sweden; youen.grusson@wanadoo.fr</footnote>

<footnote n="4">⁴ INSA, Université de Toulouse, 135 Avenue de Rangueil, 31077 Toulouse, France; rosetlaure@gmail.com</footnote>

<footnote n="5">⁵ Hydro-Environmental Processes Research Group, Applied Chemistry Department, Faculty of Chemistry, University of the Basque Country UPV/EHU, 2018 Donostia-San Sebastian, Basque Country, Spain; maite.meaurio@ehu.eus (M.M.); juanluis.lechuga@ehu.eus (J.L.L.-C.)</footnote>

<footnote n="6">⁶ Pyrenean Institute of Ecology, Spanish National Research Council, Av. Montañana 1005, 50059 Zaragoza, Spain; enrique.navarro@ipe.csic.es (E.N.); msevilla@ipe.csic.es (M.S.-C.)</footnote>

<footnote n="7">⁷ Applied Physics Department, Faculty of Marine and Environmental Sciences, Campus Puerto Real S/N, University of Cádiz, 11510 Puerto Real, Cádiz, Spain; juanjesus.gomiz@uca.es</footnote>

<footnote n="8">⁸ Department of Geological and Mining Engineering, University of Castilla La Mancha, 45002 Toledo, Spain; josemaria.bodoque@uclm.es</footnote>

<footnote>* Correspondence: cakir.roxelane@gmail.com (R.C.); sabine.sauvage@univ-tlse3.fr (S.S.); Tel.: +33-613-523-804 (R.C.); +33-643-983-625 (S.S.)</footnote>
# 1. Introduction

In the last few decades, the increase of the world population, the living standard improvement, the increase of consumption patterns and the expansion of irrigated agriculture have driven a rising global demand for water [1-3]. In Europe, the summer water scarcity impacts ecosystems and affects human management over the entire South-Western Europe [4]; (http://www.aguamod-sudoe.eu/en/). Climate change and anthropogenic activities increase the risk of water quantity and quality deterioration [5,6]. The European Union is currently supporting research programs helping regional and local governments to establish and implement policies addressing such environmental issues (i.e., http://www.aguamod-sudoe.eu/en/).

Models have become key tools to support decision-making [7,8]. Modeling is indeed a useful approach to sustainably manage water resources and assess the needs of different users over a given territory [4]. Output generated by agro-hydrological models offer many insights for both stakeholders and the scientific community, concerning water resources management, assessment of nitrate load and contamination sources, determination of ecological function's behavior or forecast of the potential impacts of global and local change, including environmental policies [9]. Deploying a model with an homogeneous resolution at the European level could provide an integrated approach of the territory and answer the European Union's need about the harmonization of transboundary tools and methodologies focused on water resources and global changes. To conduct such cross-boundary water analysis in a context of global changes, the model deployed needs to (i) incorporate climate, (ii) conceptualize different types of water resources such as "blue" (surface and aquifers) and "green" (evapotranspiration and soil moisture) waters (cf. [10,11]), (iii) have a detailed representation of limiting nutrients and sediment dynamics, and (iv) integrate anthropogenic pressures such as agricultural, industrial, and urban inputs. Investing modeling methodology at large scale is also important in order to understand water and nutrient cycles, as well as watershed behavior. The modeling methodology can be applied to any model as soon as the model fulfills large scale criteria and could allow the researchers to work on model uncertainties in order to improve the model conceptualization [12,13].

Fu et al. [8] and Moriasi et al. [13] assessed a large number of existing hydrological and water quality models and their calibration/validation procedures. According to their results, water quality models can be divided into two types regarding modeling based on purpose [8]: (1) Scientific exploration (understanding the sources/fate of biogeochemical cycles and contaminants) and (2) decision support (e.g., assessing the impacts of climate and management changes). The scales of these models can vary from field-scale (APSIM model [14]) to large scale (INCA model [15]) and the time step can vary from minutes, days [16,17], months, years [18,19], to even statics [14].

Over the past 15 years, the most commonly used hydrological model has been the Soil and Water Assessment Tool (SWAT) based on the Scopus environmental field database containing 3282 scientific papers using SWAT [8]. SWAT [16] is an open-source model with a large number of possible applications as seen in various studies [8,20,21]. SWAT can be applied from local to continental scale [9,22]. SWAT incorporates different modules issuing from modified stand-alone models such as QUAL2E and EPIC modules-developed by Brown and Barnwell [23] and Sharpley and Williams [24] respectively to simulate water quality and crop production. Moreover, management practices application is a strength of SWAT through the multiple options of inputs such as fertilizer, tillage practice, irrigation, and management practices. Finally, the SWAT model is also able to integrate urban and point sources inputs.

Developed by the United States Department of Agriculture (USDA), SWAT has been used around the world to simulate water quality and quantity [8,22], for large scale on finer resolution such as the U.S. territory with the "Hydrologic Unit Model for the United States" [25,26] or the whole Europe [9]. Abbaspour et al. [9] demonstrate the feasibility of applying SWAT at the continental scale in Europe with a subbasin resolution and at monthly time steps. At the continental scale, the main restrictions are the lack of harmonized observed data and the lack of modelers' knowledge on the functioning of all the watersheds to calibrate the model [9]. Calibration of hydrological model is based on the water balance within the watersheds. The water balance is the difference between inputs (precipitation) and outputs
(discharge and evapotranspiration) at the outlet of the watershed. However, even if the water balance is correctly simulated, it does not mean that the water in the different compartments (soil, shallow and deep aquifers, and the fluxes between these compartments) are corrected [11]. During the last decade, attempts have been made to calibrate the hydrological model on additional hydrological parameters than discharge in order to analyze the water distribution between the different compartments and fluxes inside the watershed [11]. The complexity of such calibration increases with the level of human activity. Human activities affect water resources [27,28] such as water used for irrigation, hydroelectric dam operations, or artificial pond feeding. A model does not initially integrate this level of complexity but the calibration can help to incorporate the human influence into the model. Hence, a dataset and tools to calibrate the model are therefore important if the model is to accurately represent water stocks and fluxes in the different compartments of the hydrological cycle. Such calibration approach has been shown to be necessary in order to enhance the management of water resources especially in the context of global change [29–31]. This study provides suggestions on the model set up choice and the calibration and validation procedures according to the level of the Index of Hydrological Alteration (IHA) [32] over a large-scale-transborder-European region. To achieve this goal, South West Europe was chosen as the study area due to the combination of climatic conditions and human stressors leading to recurrent water scarcity episodes [33].

The SWAT model was set up to build a calibrated and validated hydrological and water quality model over South West Europe at subbasin level and monthly time steps. Water quality simulation is evaluated through sediment load, reflecting particulate pollutants and nitrate loads, depicting dissolved pollutants [20]. The objective of this study is to evaluate the need to frame the calibration not only on the water balance knowledge but also in the water and agricultural management practices (irrigation, agricultural inputs). For example, irrigation over South West Europe is important, representing 58% of the water use in Spain and up to 90% in regions with intense agricultural practices such as Extremadura in South West Spain [34]. Several studies pointed out the issue of water scarcity in South-Western Europe [2,33,35]. It is thus essential to have a good representation of the volume of irrigation in the model, as irrigation volume will constrain water balance, especially in climate change and cultural practices impact studies. In the present study, we suggest an additional validation method of water balance through the comparisons of spatialized water irrigated or/and crop yield with national datasets. This additional validation might be useful for studies with limited documentation and a lack of financial support. One possible approach would be exploring the potential correlation between the calibration/validation procedure and the IHA [32]. The IHA characterizes the impact of river regulation from anthropogenic activity such as dams and irrigations on the flow regime compared to the natural regime. IHA is a well-known hydrological indicator used in environmental studies [36–38] to scale the anthropogenic disturbance level of a watershed. No study has yet attempted to link the level of complexity of a calibration procedure with the IHA.

Therefore, the partial objectives of this work are: (i) To calibrate and validate a multi-countries scale multi-watershed models based on different open-source data available at national and European scale, (ii) to provide a calibration guideline regarding the degree of human activities in the multi-watersheds and (iii) to determine the strength and weaknesses of this guideline.

## 2. Materials and Methods

### 2.1. Study Area

The study area, called SUDOE, covers 773,191 km² in South-Western Europe (Figure 1A). SUDOE is a transnational territory delimited by the European Union in order to develop scientific and institutional collaborations and to build integrative approaches regarding environmental issues. SUDOE encompasses the south west of France, Andorra, and the Iberian Peninsula (Spain, Portugal, and Gibraltar). More than 60 million inhabitants live in the SUDOE territory, including a few of the biggest European cities (Figure 1A). This territory is the main agricultural region of the European
Union with about half of its surface dedicated to this activity according to Corine Land Cover 2012 (CLC) map [39] (Figure 1B). The rest of the territory is occupied by 46% forests and semi-natural areas, 3% artificial surfaces and 1% wetlands and water bodies. A substantial part of the agricultural surface is irrigated; for instance, about 13.5% are permanently irrigated in Spain and Portugal [34]. The irrigation volume is important in SUDOE territory where water stress is more acute [2,33].

<figure><img src="image_4.png" /><figcaption>**Figure 1.** SUDOE study area characteristics (A) Topography, main cities and regions; (B) Land use and subsystems; (C) Pedology.</figcaption></figure>

The SUDOE two main climates are (i) Mediterranean over most parts of the Iberian Peninsula and the Mediterranean coast of France and (ii) temperate oceanic from the northern coast of Spain to the west coast and inland part of the south west of France [40]. The average annual precipitation measured between 1993 and 2012 over the area is 785.7 mm [41]. SUDOE territory incorporates major and minor streams directly flowing to the sea/ocean that corresponds to 78.75 km³·year⁻¹ of freshwater water resources over the period 1993–2012 [41] with an average nitrate concentration of 1.7 mg·L⁻¹ [4] over the period 1992–2012. Hydrological systems are managed differently by stakeholders in each subsystem depending on the administrative entity (Figure 1B).

Urban sewage has major impacts within the SUDOE surface water resources where the reduction of streamflow during summer is the main issue for health security [4]. Low streamflows increase pollutant concentrations in rivers by decreasing the dilution effect, affecting aquatic ecosystems and threatening the availability of good quality water.

## 2.2. *Global Methodology*

In order to conceptualize the water and nutrient balance of the entire SUDOE territory, the model has been divided into three zones characterized by a specific pedoclimatic condition: The Pyrenean, the Mediterranean, and the Oceanic zone (pedologic information in Figure 1C and the division of the pedoclimatic regions in Figure 2). The Pyrenean zone is characterized by an oceanic climate and the presence of the Pyrenean mountains; the Mediterranean zone is described by a Mediterranean climate and, the Oceanic zone has been identified in the west part of the territory with an oceanic climate.
<figure><img src="image_5.png" /><figcaption>**Figure 2.** SUDOE models ((1) Pyrenean project, (2) Mediterranean projet, (3) Oceanic Project; (A) Garonne, (B) Adour, (C) Aude, (D) Minho, (E) Duero, (F) Ebro, (G) Tajo, (H) Guadalquivir, (I) Jucar, (J) Guadiana, and (K) Segura).</figcaption></figure>

The model set up and calibration/validation procedures were identical for each project. The model was set up, calibrated, and validated on monthly time steps thanks to the collaborations of French, Portuguese, and Spanish stakeholders working in public (Water agencies, non-profit organizations) and private institutions (hydroelectric and dam management companies). These experts shared their data and gave their advice and suggestions on the integration of dams, crop management practices, irrigation volume and crop yield in the modeling calibration framework. The model's set up incorporates dam management, urban sewages, and crop management.

Two methodologies for calibration of the model has been done: (i) A conventional calibration (CC), based on recorded in-stream water quality and quantity integrating auto-fertilization and auto-irrigation, i.e., the model estimates the amount of fertilizer and irrigation demanded by plants, which is the most common calibration strategy for water resource modeling and (ii) an additional calibration (AC), based on the same data, but including additional outputs coming from other models' outputs or from national statistics (e.g., irrigation and crop yield, Figure 3). The comparison of these two calibrations allows this study to determine the strengths and weaknesses of the CC and AC methodologies on the representation of water resources. The last step of this study is to define an index that indicates the choice of the calibration/validation procedure.

## 2.3. Modeling Approach: The SWAT Model

In the SWAT model, the simulated watersheds are first divided into subbasins based on river network and topography. The implementation procedure allows us to identify Hydrological Response Units (HRUs) in each subbasin based on a unique combination of soil, land cover and slope information. HRUs are used to compute water balance articulated around four compartments (soil, snow, shallow, and deep aquifer) linked by hydrological processes such as infiltration, runoff, evapotranspiration, lateral flow, and percolation. HRUs water balance is then summed at the subbasin scale and routed to the watershed outlet. Impact on hydrology, sediment and nitrate dynamics of urban sewage, agricultural management and dams are taken into consideration at subbasin scale.
Modeling SUDOE territory

<figure><img src="image_6.png" /><figcaption>**Figure 3.** Global approach of this study to model SUDOE territory.</figcaption></figure>

The surface runoff volume is computed from a modified version of the Soil Conservation Service (SCS) curve number procedure [42]. Muskingum routing method was selected to calibrate streamflow in the channel [43]. The Hargreaves method [44] was chosen to estimate the potential evapotranspiration, as the Hargreaves method only requires air temperature and precipitation. Hargreaves's formula has been shown to satisfactorily simulate the potential evapotranspiration when applied outside of humid climatic zones [45].

In the SWAT model, sediment transport is simulated both in-land and channel compartments. Over lands, erosion and sediment yield is estimated using the Modified Universal Soil Loss Equation (MUSLE) [46] for each HRU. This equation considers hydrology and land characteristics (more details of the USLE equation factors can be found in Neitsch et al. [47]).

Inland and instream nitrogen cycles are also modeled. Soil and atmospheric waters, as well as soil characteristics, drive the nitrogen cycle within each HRU. Nitrogen can be added to the system by fertilizer applications, dry deposition, and rainfall. Instream nitrate processes are based on the QUAL2E model integrated into the SWAT model [23,48]. More details of the SWAT nitrate module can be found in Brown and Barnwell [23] and Neitsch et al. [44].

The plant growth component of SWAT is based on the radiation use efficiency approach with empirical parameters. Plant growth can be inhibited by temperature, water, nitrogen and phosphorus stress. Plant development is based on daily accumulated heat unit values [49,50]. The heat unit states the stage of plant development. It varies from 0 to 1, 0 indicating the sowing time and 1 the optimal moment for the plant to be harvested.

Potential biomass is estimated using a method developed by Monteith et al. [51] depending on leaf development, light interception, and light-biomass ratio conversion. This potential is calculated considering ideal growing conditions with adequate water and nutrient supply and a favorable climate.

When a harvest/kill operation occurs, a portion of plant biomass is removed from the system as yield. Crop yield is calculated from biomass production, harvest efficiency and the fraction of the above-ground plant's dry biomass removed as dry yield also called harvest index. Harvest index is calculated daily according to plant stage development and HRU water deficit.
## 2.4. Input Data

A large database compilation has been gathered to study the transnational SUDOED area. Data sources, as well as their spatial and temporal resolutions, are displayed in Table 1. CC and AC are using the same model set up and have the same input data. A global 1-km resolution land surface Digital Elevation Model (DEM), derived from the United States Geological Survey (USGS), was used to described topography and quantify slopes. Such DEM has been selected as a compromise between the large-scale model setup and the computation resources available. Land cover was available from Corine Land Cover map 2012 [39]. Soil characteristics were defined from the Harmonized world soil map that encompasses the dataset of the European Soil Database.

<table><caption>**Table 1.** Data sources.</caption><thead><tr><th rowspan="2">Data Type</th><th rowspan="2">Data Source</th><th colspan="2">Resolution</th></tr><tr><th>Spatial</th><th>Temporal</th></tr></thead><tbody><tr><td>DEM</td><td>Shuttle Radar Topography Mission [52]</td><td>~1 km</td><td>Static based in 2006</td></tr><tr><td>Land cover</td><td>Corine Land Cover 2012 [39]</td><td>100 m</td><td>Static based in 2012</td></tr><tr><td>Soil</td><td>Harmonized world soil database [53,54]</td><td>~1 km</td><td>Static based in 2012</td></tr><tr><td>Climate (Observed)</td><td>MESAN reanalyses of precipitation and temperature (SMHI) [55,56]</td><td>~12 km</td><td>Daily (2000–2010)</td></tr><tr><td>City effluents</td><td>UWWTP—EUDB (https://ec.europa.eu/)</td><td>Vectorial points</td><td>Annual (2010)</td></tr><tr><td rowspan="3">Crop management</td><td>France: AGRESTE (http://agreste.agriculture.gouv.fr/)</td><td>Regional</td><td>Annual (2017)</td></tr><tr><td>Spain: MAPAMA (https://www.mapama.gob.es/)</td><td>Regional</td><td>Annual (2015)</td></tr><tr><td>Portugal: INE (https://ine/pt)</td><td>Regional</td><td>Annual (1998)</td></tr><tr><td rowspan="2">Irrigation</td><td>France: CACG (https://www.cacg.fr/fr/)</td><td>Regional</td><td>Monthly (2016)</td></tr><tr><td>Spain: AQUATOOL (https://aquatool.webs.upv.es/)</td><td>Subsystems</td><td>Monthly (2000–2010)</td></tr><tr><td rowspan="2">Dam management</td><td>France: EDF (REGARD-RTRA/STAE)</td><td>Vectorial points</td><td>Annual (1923–2017)</td></tr><tr><td>Spain: MAPAMA (https://www.mapama.gob.es/)</td><td>Vectorial points</td><td>Annual (1944–2017)</td></tr><tr><td rowspan="5">River discharge</td><td>France: Banque Hydro (http://www.hydro.eaufrance.fr/)</td><td>Vectorial points</td><td>Daily (1921–2013)</td></tr><tr><td>Spain: MAPAMA (https://www.mapama.gob.es/)</td><td>Vectorial points</td><td>Daily (1911–2017)</td></tr><tr><td>Spain: SIMPA “natural” discharge simulation (http://www.cedex.es/)</td><td>Vectorial point</td><td>Daily (1980–2010)</td></tr><tr><td>Portugal: SNIRH (https://snirh.apambiente.pt/)</td><td>Vectorial points</td><td>Daily (2000–2014)</td></tr><tr><td>France: Agence de l'eau (http://adour-garonne.eaufrance.fr/)</td><td>Vectorial points</td><td>Punctual (1980–2016)</td></tr><tr><td rowspan="2">Nitrate</td><td>Spain: Confederaciones (https://http://hispagua.cedex.es/instituciones/confederaciones/)</td><td>Vectorial points</td><td>Punctual (1971–2014)</td></tr><tr><td>France: Agence de l'eau (http://adour-garonne.eaufrance.fr/)</td><td>Vectorial points</td><td>Punctual (1980–2016)</td></tr><tr><td rowspan="2">Sediment</td><td>Spain: Confederaciones (https://http://hispagua.cedex.es/instituciones/confederaciones/)</td><td>Vectorial points</td><td>Punctual (1971–2014)</td></tr></tbody></table>

Daily precipitation and temperature originate from the MESoscale ANalysis (MESAN) climate reanalyses spanning Europe for the 1989–2010 period [55,56]. This dataset has been shown to be as good as the reference French meteorological dataset [57], SAFRAN [58]. In Spain, the MESAN reanalysis is not based on a dense network of observations, but a comparison of streamflow outputs forced by SAFRAN [59] and MESAN in different rivers over Spain showed similar results.

### 2.4.1. Dam Management

The characteristics of dams (volume, year of construction, water release, management rules) were collected for Spain and France (Table 1). Only dams with a total volume superior to 5 Hm³ were
considered to strongly impact the hydrology and taken into consideration. In the case of different dams in the same subbasin, a fictional dam conceptualization integrating the volume equivalent of the sum of dams has been used. A total of 44 dams were implemented (Figure 2). Management data were used to set minimum and maximum objective volumes or the target storage water volume. In some cases, the maximum and minimum target water release was reconstructed from the observed streamflow downstream of the reservoir.

### 2.4.2. Urban Sewage

Urban sewage localizations and characteristics (maximum capacity, annual load) were taken from the 2010 European dataset [60] including information about all European cities. Every wastewater treatment plant (WWTP) with more than 50,000 equivalent inhabitants' capacity has been considered. Such as for dams, the WWTPs were aggregated by subbasin. A total of 44 point sources were integrated into the model (Figure 2). The annual spilled volume is deduced from the number of equivalent-inhabitants of selected WWTPs and the annual water consumption per inhabitant per country (France: 58.04 m³ year⁻¹, Portugal: 67.46 m³ year⁻¹ and Spain: 47.4 m³ year⁻¹; https://www.insee.fr, https://ine/pt, https://ine/es). In order to overcome the lack of nutrient loading from WWTPs at the European scale, a simplified version of the Zessner and Lindter [61] method has been used to estimate the annual nitrate loading with the available information about treatment plant. The daily averaged nitrogen loads that enter the treatment plant is assumed to be 8.8 gN person-equivalent⁻¹ and the average treatment efficiency for the whole South West of Europe is 40% [61]. The following equation determines the load of nitrogen released in the environment:

$$X_{load} = (1 - X_{removal}) \times Load_{design} \times X_{influent load} \quad (1)$$

with:

$X_{load}$: The nitrogen pollution in g,

$X_{removal}$: The fraction of efficiency between 0 and 1,

$Load_{design}$: The inhabitant equivalent,

$X_{influent load}$: The daily average amount of nitrogen rejected per inhabitant (g/inhabitant equivalent/d).

This equation has been validated at the outlet of several French and Spanish watersheds, showing a good correlation between observed and simulated nitrogen loads (R² = 0.87, p-value < 0.01) (Supplementary Materials Figure S1).

### 2.4.3. Crop Management

Interannual average of monthly irrigation volumes for the South of France and Spain were given by the Institut national de recherche en sciences et technologies pour l'environnement et l'agriculture (IRSTEA) [62–64], the Compagnie d'Aménagement des Coteaux de Gascogne (CACG) and the Spanish Hydrological National Plan [65,66]. Crop management data were obtained from the national database of agriculture statistics from each country. The technical management and fertilizer application for each crop in each region originate from the recommendations provided by the French Chambers of Agriculture (FCA), the French and Spanish Directorate-General for Agriculture and Rural Development (DGADR) and the Spanish Confederations (SC).

Crop management were applied to HRUs according to the recommendations of FCA, DGADR, and SC. The dominant crops (determine by national statistics) of each region and crop technical management were selected and implemented in the model used for both CC and AC calibrations (Table 2). For example, the dominant irrigated crop is cotton in the South East of Spain whereas in France dominant irrigated crops are corn and wheat.
<table><caption>**Table 2.** Type of crops and technical management integrated into the model for CC and AC calibrations recommended by institutions. For AC, the amount of fertilizers was calibrated according to the crop yields outputs (Supplementary Materials Table S1).</caption><thead><tr><th>Territory</th><th>Country</th><th>Classification Description from CLC 2012</th><th>CLC Code</th><th>Selected Crops</th><th>Fertilizers Type</th><th>Quantity Applied (kg ha<sup>-1</sup> year<sup>-1</sup>)</th></tr></thead><tbody><tr><td rowspan="10">France</td><td rowspan="10">France</td><td>Permanently irrigated land</td><td>2.21.212</td><td>Maize</td><td>15-15-15 <sup>(1)</sup></td><td>400</td></tr><tr><td>Non-irrigated arable land</td><td>2.21.211</td><td>Wheat</td><td>15-30-0 <sup>(2)</sup></td><td>200-150</td></tr><tr><td>Land principally occupied by agriculture, with significant areas of natural vegetation</td><td>2.24.243</td><td>Maize</td><td>15-15-15 <sup>(1)</sup></td><td>400</td></tr><tr><td>Vineyard</td><td>2.22.221</td><td>Vine</td><td>N<br>P</td><td>30<br>20</td></tr><tr><td>Fruits trees and berry plantation</td><td>2.22.222</td><td></td><td>18-46-0 <sup>(3)</sup></td><td>100</td></tr><tr><td>Rice fields</td><td>2.21.213</td><td>Rice</td><td>N<br>P<br>K</td><td>50<br>60-80<br>60-100</td></tr><tr><td>Non-irrigated arable land</td><td>2.21.211</td><td>Wheat</td><td>N<br>P<br>K</td><td>240<br>96<br>144</td></tr><tr><td>Land principally occupied by agriculture, with significant areas of natural vegetation</td><td>2.24.243</td><td>Sunflower</td><td>15-15-15 <sup>(1)</sup></td><td>400-600</td></tr><tr><td>Permanently irrigated land</td><td>2.21.212</td><td>Cotton</td><td>N<br>K</td><td>150<br>33</td></tr><tr><td rowspan="10">Pyrenean</td><td rowspan="10">Spain, Portugal</td><td>Vineyard</td><td>2.22.221</td><td>Grape</td><td>N<br>P</td><td>30<br>20</td></tr><tr><td>Fruits trees and berry plantation</td><td>2.22.222</td><td></td><td>N<br>P<br>K</td><td>190<br>120<br>360</td></tr><tr><td>Rice fields</td><td>2.21.213</td><td>Rice</td><td>N<br>P<br>K</td><td>190<br>60<br>240</td></tr><tr><td>Almonds</td><td>2.22.222</td><td>Almonds</td><td>N<br>P<br>K</td><td>245<br>112<br>156</td></tr></tbody></table>
<table><caption>Table 2. *Cont.*</caption><thead><tr><th>Territory</th><th>Country</th><th>Classification Description from CLC 2012</th><th>CLC Code</th><th>Selected Crops</th><th>Fertilizers Type</th><th>Quantity Applied (kg ha<sup>-1</sup> year<sup>-1</sup>)</th></tr></thead><tbody><tr><td rowspan="12">Oceanic</td><td rowspan="12">Spain, Portugal</td><td rowspan="3">Non-irrigated arable land</td><td rowspan="3">2.21.211</td><td rowspan="3">Wheat</td><td>N</td><td>290</td></tr><tr><td>P</td><td>126</td></tr><tr><td>K</td><td>144</td></tr><tr><td rowspan="2">Land principally occupied by agriculture, with significant areas of natural vegetation</td><td rowspan="2">2.24.243</td><td rowspan="2">Cotton</td><td>N</td><td>150</td></tr><tr><td>P</td><td>33</td></tr><tr><td rowspan="2">Permanently irrigated land</td><td rowspan="2">2.21.212</td><td rowspan="2">Cotton</td><td>N</td><td>150</td></tr><tr><td>P</td><td>33</td></tr><tr><td rowspan="2">Vineyard</td><td rowspan="2">2.22.221</td><td rowspan="2">Grape</td><td>N</td><td>60</td></tr><tr><td>P</td><td>40</td></tr><tr><td rowspan="3">Olive groves</td><td rowspan="3">2.22.223</td><td rowspan="3">Olive</td><td>N</td><td>133</td></tr><tr><td>P</td><td>52</td></tr><tr><td>K</td><td>195</td></tr><tr><td rowspan="3">Rice fields</td><td rowspan="3">2.21.213</td><td rowspan="3">Rice</td><td>N</td><td>190</td></tr><tr><td>P</td><td>60</td></tr><tr><td>K</td><td>240</td></tr><tr><td rowspan="3">Almonds</td><td rowspan="3">2.22.222</td><td rowspan="3">Almonds</td><td>N</td><td>245</td></tr><tr><td>P</td><td>112</td></tr><tr><td>K</td><td>156</td></tr></tbody></table>
<table><caption>Table 2. *Cont.*</caption><thead><tr><th>Territory</th><th>Country</th><th>Classification Description from CLC 2012</th><th>CLC Code</th><th>Selected Crops</th><th>Fertilizers Type</th><th>Quantity Applied (kg ha<sup>-1</sup> year<sup>-1</sup>)</th></tr></thead><tbody><tr><td rowspan="12">Mediteranean</td><td rowspan="12">Spain, Portugal</td><td rowspan="3">Non-irrigated arable land</td><td rowspan="3">2.21.211</td><td rowspan="3">Melon</td><td>N</td><td>85</td></tr><tr><td>P</td><td>45</td></tr><tr><td>K</td><td>95</td></tr><tr><td rowspan="2">Land principally occupied by agriculture, with significant areas of natural vegetation</td><td rowspan="2">2.24.243</td><td rowspan="2">Lemon</td><td rowspan="2">21-0-0+24S <sup>(4)</sup></td><td rowspan="2">610</td></tr><tr><td>245</td></tr><tr><td rowspan="2">Permanently irrigated land</td><td rowspan="2">2.21.212</td><td rowspan="2">Almonds</td><td>N</td><td>112</td></tr><tr><td>P</td><td>156</td></tr><tr><td rowspan="2">Vineyard</td><td rowspan="2">2.22.221</td><td rowspan="2">Grape</td><td>N</td><td>60</td></tr><tr><td>P</td><td>40</td></tr><tr><td rowspan="3">Olive groves</td><td rowspan="3">2.22.223</td><td rowspan="3">Olive</td><td>N</td><td>133</td></tr><tr><td>P</td><td>52</td></tr><tr><td>K</td><td>195</td></tr><tr><td rowspan="3">Rice fields</td><td rowspan="3">2.21.213</td><td rowspan="3">Rice</td><td>N</td><td>190</td></tr><tr><td>P</td><td>60</td></tr><tr><td>K</td><td>240</td></tr><tr><td rowspan="3">Fruits trees and berry plantation</td><td rowspan="3">2.22.222</td><td rowspan="3">Almonds</td><td>N</td><td>210</td></tr><tr><td>P</td><td>49</td></tr><tr><td>K</td><td>190</td></tr></tbody></table>

<footnote>(1) 15-15-15, NPK fertilize, is a combination of nitrogen (N), phosphorus (P) and potassium (K) designed to maximize crop yield and quality. It contains 15%N, 15% P, and 15% K. (2) 15-30-0 is a high phosphate fertilizer. It contains 15%N, 30% P and 0% K. (3) 18-46-0, Diammonium phosphate, is a fertilizer and when applied as food for plants, the soil pH is increased temporarily. It contains 18%N, 46% P₂O₅, and 0% K. (4) 21-0-0+24S, the Ammonium Sulfate Fertilizer contains 21%N and 24% of sulfur.</footnote>
The same crop management has been applied to all crops in the SUDOE territory. The beginning of plant growth and the auto-irrigation module are set to 0.15 heat units. Auto-irrigation activation terms are defined according to the type of crop and the climatic condition. The auto-irrigation can be triggered by two conditions, chosen for each subbasin regarding calibration, based on (1) plant water demand or (2) soil water content. Due to the lack of precise data, choices on the repartition of irrigated and non-irrigated crops have been made: Rice, olive, vineyard, fruit trees, berry plantations, maize, and almonds are considered as irrigated crops, while the rest are considered as non-irrigated. Irrigation is supplied by the nearest dams or main rivers. Fertilizers are applied automatically according to plant needs. The harvest is set at 1.2 heat units. The harvest heat unit is superior to 1 because this study considers that there is always a time lag between the optimized harvest time and the true harvest [47]. In that case, between 1 and 1.2 heat units, the plant consumes energy without producing more.

The model used to calibrate CC and AC has the same setup. For both calibrations, the irrigation was set up as "auto-irrigation". For CC and AC calibration we enter a volume at each auto-application and a maximum volume. The exact volume to apply is determined by the model according to the parametrization (plant water demand or soil water content). However, for the AC, crop management routing, simulated irrigation volume, crop yield, and biomass were calibrated and validated with national statistic datasets (Table 1) whereas, for CC, these model outputs were not calibrated. By default, SWAT yield estimates are based on average yields of US crop production. Some parameters of plant growth were adjusted to South-Western European crops based on Institute of Agriculture recommendations (cf. Supplementary Materials Table S1).

## 2.5. Study Area Discretization

The conceptualized approach is represented in Figure 3. As explained in Section 2.2 Global Methodology, three projects have been set up to be able to depict three contrasting territories: The Pyrenean, the Mediterranean and the Oceanic zones (Figure 2). The delineation, warm-up and calibration/validation steps were achieved following the same procedure for each project. HRUs have been defined using a total of 24 land classes [39], 34 soil classes [67] and 5 slopes classes (0–2, 2–5, 5–15, 15–25, and >25°) [68]. A total of 263 subbasins with 22,997 HRUs have been identified in the SUDOE territory.

## 2.6. Calibration, Multiple Validations and Uncertainties

### 2.6.1. The Conventional Calibration (CC)

For the CC calibration, hydrology, sediments and nitrate loads were calibrated. A three year warm-up period from 1997 to 1999, was taken into consideration to set up the initial variables' value. Hydrology and nitrate loads were automatically calibrated based on sensitivity analysis at 47 gauging stations including 19 hydrological stations and 28 physical-chemical stations over the period 2000–2005. *p*-value (determines the significance of the sensitivity) and t-stat (provides a measure of sensitivity) were calculated to identify the sensitive parameters [69]. The LH-OAT sensitivity analysis tool [70] was used in order to determine the sensitive parameters that are presented in Table S1. The SUFI-2 algorithm of SWAT-CUP [13,69,71] was selected to perform the automatic multi-site calibration. The SUFI-2 algorithm is known to identify an appropriate parameter set in a limited number of iterations [72]. A total of 1500 runs have been done for calibration steps as recommended by Yang et al. [73]. Validation was then conducted for the same stations over the period 2006–2010. The model was run on a monthly time step.

Four standard statistical criteria suggested by Moriasi et al. [13,74] and Gupta et al. [75] were used to evaluate the goodness of hydrological simulations: coefficient of determination (R²), Nash and Sutcliff Efficiency (NSE), percent bias (PBIAS) and Kling-Gupta efficiency (KGE). The NSE is commonly used to evaluate hydrological models. It ranges from minus infinity to 1 (perfect fit). The correlation coefficient (R²) measures the fit between observed and simulated values and varies from
0 to 1 (strong linear association). A perfect fit requires that the regression slope and intercept are equal to 1 and 0 respectively. The PBIAS is a metric indicating and global underestimation (positive value) or an overestimation (negative value) of the simulation compared to observations. According to Moriasi et al. [13], NSE, R², and PBIAS are satisfactory for values higher than 0.5, 0.6, and ±15% respectively for monthly time series. KGE is a hydrological criteria suggested by Gupta et al. [75] and modified by Kling et al. [76], based on “r” Pearson coefficient (similar to R²), model bias and variability. It ranges from minus infinity to 1 (strong linear association). According to Knoben et al. [77], simulation between -0.41 and 1 could be considered as reasonable performance.

The nitrate and sediment load simulations were compared to sampled data collected once or twice a month from 2000 to 2010 (Table 1). The Load Estimator model (LOADEST) [78] was used for each physicochemical station due to the lack of nitrate and sediment data in order to generate monthly nitrate and sediment load and extend time series, as commonly used to calibrate and validate SWAT outputs [79–81]. This software generates water quality parameters through regressions between streamflow and water quality [78]. A statistically-based optimization procedure identifies the most accurate of the eleven regressions. The Spearman's rank correlation coefficient (ρ, [82]), the non-parametric test measuring the dependence between two variables, is used to compare LOADEST simulation with observations. LOADEST simulations give sediment loads (ρ = 0.94 and *p*-value < 0.001) and nitrate loads (ρ = 0.99 and *p*-value < 0.001) close to observations (Figures S2 and S3 in Supplementary Materials), even if the LOADEST model slightly underestimates the sediment load (slope = 0.71). LOADEST provides us interpolated water quality data that were used to calibrate and validate the SWAT model on a monthly scale in addition to the observed data. The R² statistic is used to evaluate the agreement of the nitrate and sediment load simulation between LOADEST and SWAT models for calibration (2000–2005) and validation periods (2006–2010).

### 2.6.2. The Additional Calibration (AC)

From the CC procedure, an additional calibration (AC) we made by adding two calibrations steps: (1) the calibration of the “natural streamflow”—the streamflow that should occur without any anthropogenic pressures and (2) the calibration of the anthropogenic pressures (Fertilization application, crop yield, and irrigation volume) in order to improve the evapotranspiration and the uptake of plants. The validation step is the same procedure as CC using SWAT-CUP. Data used to validate AC were the “natural streamflow”, the crop management data and the one used during CC procedure. As a reminder, hydrology, sediments and nitrate loads were calibrated (2000–2005) and validated (2006–2010).

To ensure the correct hydrological performance of the SWAT projects (Figure 2) calibrated with AC methodology, discharge outputs were compared with two different sources in annual and monthly time steps. The first comparison was done with the observed discharge data at the subsystem scale (administrative entities, Figure 1B). The second one was performed with the discharge values simulated with the “Integrated System for Rainfall-Runoff Modeling” model, commonly named SIMPA [83]. The AC “natural” streamflow simulated by SWAT (the simulation where all human pressures were withdrawn) was calibrated by comparing it with SIMPA simulations, which simulated natural flows. SIMPA is a rainfall-runoff model, developed by the Center for Hydrographic Studies of the Spanish Ministry of Environment (CEDEX) to analyze the spatial and temporal hydrological variables and to simulate hydrological processes without anthropogenic pressures. In Spain, SIMPA model is a reference for public and private stakeholders [84,85]. To compare our SWAT model with SIMPA, the CC and AC calibrations were run after deleting all anthropogenic pressures (agriculture, urban inputs and dam management). Inter-annual and monthly “natural” water volumes were compared with SIMPA outputs in the 7 SUDOE main watersheds (Figure 2). At each calibration step, hydrology, sediments, nitrate loads, and dam and crop managements were re-verified and re-calibrated if necessary.

Dam management, fertilization application, crop yield and irrigation volume were calibrated from 2000 to 2005 and validated from 2006 to 2010. Simulated irrigation volumes were compared to
the annual national values. In Spain, where data were available, a comparison of irrigated volume at the regional scale has also been conducted. Moreover, a comparison between simulated crop yields, fertilizer amount applied and national statistics was made on an annual basis. To quantify the goodness of fit between simulated and observed mean crop yield, a PBIAS value is calculated, defined by:

$$ PBIAS = \frac{\sum_{1}^{n}(SWAT simulated yield - National yield) * 100}{\sum_{1}^{n} National yield} \quad (2) $$

where n is the number of observations in the period under consideration.

## 2.7. Climate Change and Agricultural Practices Changes

In order to investigate the sensitivity of the model to climate change, the variability of the simulated outputs has been tested for a precipitation decrease of 20% over the period 1980–2010 based on the global precipitation change between 1971 and 2000 and the +2 °C period of Vautard et al. study [86]. CC and AC calibrations with observed precipitation have been compared to outputs where the precipitations have been decreased between 2000 and 2010. Focus has been placed on water yield (mm), fertilization (kg ha⁻¹), biomass (t ha⁻¹), yield (t ha⁻¹), and irrigation volume (Hm³·year⁻¹) and a percentage of variation has been calculated for the outputs simulated by the reduced precipitation with respect to the normal precipitation dataset. The water yield is considered as the long-term average flow to convey the idea that only a volume of water is being referred to, as opposed to a hydrograph [87].

## 2.8. Indicators of Hydrological Alteration

More than 170 hydrological indicators have been developed in the past years including the widely used Index of Hydrological Alteration (IHA) [36]. IHA characterizes the impact of river regulation on flow regimes in environmental flow studies [36]. The IHA was developed by the Nature Conservancy organization [32,37,38]. The IHA calculated on a monthly scale is the sum of twelve indicators that reflects, from different statistics (magnitude index, variability index, flood and drought intensity and frequency indexes), the modification of the current human-impacted flow regime, with respect to the natural one [32,37,38]. The IHA used in this study was proposed by Richter et al. [32] and modified by Santa-María [88]. The IHA varies between 0 and 100 with the maximum value meaning that the river flow is not altered.

In order to provide suggestions on the model set up choice and the calibration/validation procedures according to the level of the IHA, an inter-model comparison between AC and CC has been conducted. The IHA has been compared to the deviation of the hydrological criteria KGE (ΔKGE) between the CC and AC simulations. KGE has been chosen in this analysis because KGE simultaneously encompasses the model correlation, bias and variability. The ΔKGE were calculated for the 47 gauging stations considered in this study (even if the streamflows have been calibrated on only 19 gauging stations). When the ΔKGE obtain from the comparison of the CC and AC calibrations is close to 0, AC and CC performances (the bias, dispersion and correlations compared to observations) are similar.

# 3. Results and Discussions

Calibration is a crucial step on modeling, considering its impacts on model outputs and result interpretations [89]. For hydrological models, only the streamflow and the water balance is usually calibrated based on field records and previous studies [11,20,80,90]. For agronomic models, biomass production and crop yield help to calibrate the model [91]. In this study, the main difficulty regarding the modeling set up and calibration methodology is related to the large scale. Collecting data (biomass production, plant growth stage, water and nitrate runoff) is relatively easy to conduct on a local scale but becomes very complicated with time and budget constraints when a transnational scale is
considered. The results and implications of the calibration strategy for water quantity and quality are
detailed and discussed below.

3.1. *Water Quantity Evaluation*

To evaluate water quantity, the comparison of streamflow recorded in gauging stations and streamflow simulations is commonly done. This section presents the streamflow evaluation of CC and AC simulations and a comparison of these simulations. Moreover, a comparison of the water yield was conducted in order to validate the hydrological models not only on the streamflow record.

3.1.1. Streamflow Evaluation of Conventional Calibration (CC)

The simulated monthly streamflow from the conventional calibration (CC) (calibrated only on streamflows, nitrate, and sediment loads) is shown in Figure 4. The sensitivity analysis for streamflow identified 15 most influential parameters for the 45 selected (cf. Supplementary Materials Tables S2 and S3). The sensitive parameters ranking change between the different parts of the SUDOE territory. The water quantity modeling is influenced the most by the base flow alpha-factor and the groundwater delay. According to Moriasi et al. [13,74], NSE are poor (−0.3) to good (0.71), R² and PBIAS are poor to good (respectively (0.58; 0.95) and (−47.06; 4.51) and KGE are satisfactory (0.01; 0.73) for calibration period according to Knoben et al. [77]. For validation period, the range of NSE, R², PBIAS and KGE are respectively (−0.48; 0.80), (0.39; 0.95), (−41.51; 3.83), and (0.1; 0.79) (Table 3). The results are not satisfactory in terms of extreme values, however, the trends were adjusted appropriately. The global behavior of the multiple watersheds model is satisfactory to good, especially during low flow periods (Figure 4) considering the complexity of a multi-watershed model. NSE indicator is poor in some watersheds that are subjected to extreme values (e.g., Duero in 2010, Guadiana in 2001 and 2010) or highly anthropized (e.g., Jucar) (Figure 4). There are two major reasons for the poor NSE performances in some watershed: (1) reanalysis of MESAN as input climate datasets and (2) level of human stressors. NSE performances at Adour, Duero, Jucar, and Minho could be related to high uncertainties on the MESAN reanalyzes over the Pyrenees and over the South-Western coasts of Europe such as Aquitaine, Basque Country, Tajo, Duero, and Jucar regions [56].

In anthropized watersheds, the SWAT model has the ability to integrate human stressors [21] even if some management practices cannot be reproduced by the model. For example, dam management in the Tajo river watershed significantly changes between the calibration and the validation period. Streamflow does not go above 40 m³/s after 2005 due to dam management, where the peak flows are erased, whereas SWAT simulations simulate these peak flows. In addition, the multiple watershed modeling set up presented in this study has required several simplifications leading to an increase in the uncertainties. For instance, the dominant crop of each region and its technical management were selected in order to simplify the model conceptualization. If we compare this study simulation with the simulation of Abbaspour et al. [9], simulations are in the same range or better for our simulation. However, the South of Spain has poor simulation in Abbaspour et al. [9] in Jucar and Guadalquivir outlets (R² < 0.1) as much as our study (R² equal to 0.58 for Jucar outlet and to 0.67 for Guadalquivir outlet for the calibration period but NSE equal to −0.3 and 0.27 for Jucar and Guadalquivir respectively).
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Streamflow time series between simulations, conventional calibration (CC, dotted lines) and additional calibration (AC, straight lines) and observations (points) at monthly scale at the main outlets for calibration (2000–2005) and validation (2006–2010) periods (cf. Figure 2 for watershed localization: (A)—Garonne, (B)—Adour, (C)—Aude, (D)—Minho, (E)—Duero, (F)—Ebro, (G)—Tajo, (H)—Guadalquivir, (I)—Jucar, (J)—Guadiana, and (K)—Segura).</figcaption></figure>
<table><caption>**Table 3.** Nash and Sutcliff Efficiency (NSE), R², percent bias (PBIAS), and Kling-Gupta efficiency (KGE) performance of Soil and Water Assessment Tool (SWAT) models at watersheds' outlets (cf. Figure 2 for station localization) for calibration (2000–2005) and validation period (2006–2010) on a monthly scale. Bold values indicate the best performances (CC: conventional calibration; AC: additional calibration).</caption><thead><tr><th colspan="2">Models</th><th colspan="6">Conventional Calibration (CC)</th><th colspan="6">Additional Calibration (AC)</th></tr><tr><th rowspan="3">Sigle</th><th rowspan="3">Watershed</th><th colspan="4">Calibration</th><th colspan="4">Validation</th><th colspan="4">Calibration</th><th colspan="4">Validation</th></tr><tr><th colspan="4">2000–2005</th><th colspan="4">2006–2010</th><th colspan="4">2000–2005</th><th colspan="4">2006–2010</th></tr><tr><th>NSE</th><th>R<sup>2</sup></th><th>PBIAS</th><th>KGE</th><th>NSE</th><th>R<sup>2</sup></th><th>PBIAS</th><th>KGE</th><th>NSE</th><th>R<sup>2</sup></th><th>PBIAS</th><th>KGE</th><th>NSE</th><th>R<sup>2</sup></th><th>PBIAS</th><th>KGE</th></tr></thead><tbody><tr><td>B</td><td>Adour</td><td>0.10</td><td>0.92</td><td>-12.73</td><td>0.28</td><td>0.08</td><td>0.91</td><td>-14.36</td><td>0.39</td><td>0.44</td><td>0.81</td><td>10.36</td><td>0.70</td><td>0.49</td><td>0.77</td><td>-28.18</td><td>0.63</td></tr><tr><td>C'</td><td>Agly</td><td>0.24</td><td>0.77</td><td>-18.03</td><td>0.45</td><td>-0.18</td><td>0.82</td><td>-41.51</td><td>0.31</td><td>0.51</td><td>0.79</td><td>-10.36</td><td>0.65</td><td>0.49</td><td>0.87</td><td>-15.29</td><td>0.61</td></tr><tr><td>C</td><td>Aude</td><td>0.39</td><td>0.94</td><td>-2.77</td><td>0.61</td><td>0.38</td><td>0.90</td><td>-5.42</td><td>0.46</td><td>0.72</td><td>0.88</td><td><strong>2.38</strong></td><td>0.37</td><td>0.73</td><td>0.91</td><td>-11.93</td><td>0.82</td></tr><tr><td>E</td><td>Duero</td><td>0.22</td><td>0.86</td><td>-4.40</td><td>0.53</td><td>0.40</td><td>0.67</td><td>3.83</td><td>0.36</td><td>0.37</td><td>0.85</td><td>-3.97</td><td>0.56</td><td>0.36</td><td>0.68</td><td>4.30</td><td>0.33</td></tr><tr><td>F</td><td>Ebro</td><td>0.39</td><td>0.90</td><td>-4.47</td><td>0.48</td><td>0.45</td><td>0.88</td><td>-6.20</td><td>0.32</td><td>0.54</td><td>0.87</td><td>-11.96</td><td>0.69</td><td>0.43</td><td>0.80</td><td>-4.10</td><td>0.53</td></tr><tr><td>A</td><td>Garonne</td><td><strong>0.71</strong></td><td>0.90</td><td><strong>-0.43</strong></td><td><strong>0.73</strong></td><td><strong>0.80</strong></td><td>0.94</td><td><strong>-3.65</strong></td><td><strong>0.79</strong></td><td><strong>0.65</strong></td><td>0.88</td><td>3.77</td><td><strong>0.73</strong></td><td><strong>0.78</strong></td><td>0.93</td><td><strong>0.58</strong></td><td><strong>0.81</strong></td></tr><tr><td>H</td><td>Guadalquivir</td><td>0.27</td><td>0.67</td><td>-26.44</td><td>0.10</td><td>0.60</td><td>0.65</td><td>-30.39</td><td>0.35</td><td>0.58</td><td>0.53</td><td>-7.18</td><td>0.28</td><td>0.28</td><td>0.61</td><td>-5.95</td><td>0.27</td></tr><tr><td>J</td><td>Guadiana</td><td>0.35</td><td>0.67</td><td>4.51</td><td>0.14</td><td>0.17</td><td>0.39</td><td>-4.63</td><td>0.25</td><td>0.34</td><td>0.69</td><td>5.11</td><td>0.15</td><td>0.09</td><td>0.34</td><td>6.00</td><td>0.33</td></tr><tr><td>I</td><td>Jucar</td><td>-0.30</td><td>0.58</td><td>-35.80</td><td>0.34</td><td>-0.48</td><td>0.63</td><td>-32.06</td><td>0.38</td><td>0.22</td><td>0.56</td><td>-29.47</td><td>0.43</td><td>-0.06</td><td>0.61</td><td>-16.36</td><td>0.33</td></tr><tr><td>D</td><td>Minho</td><td>0.11</td><td><strong>0.95</strong></td><td>-47.06</td><td>0.52</td><td>0.07</td><td><strong>0.95</strong></td><td>-35.28</td><td>0.62</td><td>0.55</td><td><strong>0.94</strong></td><td>-38.91</td><td>0.59</td><td>0.47</td><td><strong>0.94</strong></td><td>-28.21</td><td>0.64</td></tr><tr><td>G</td><td>Tajo</td><td>0.26</td><td>0.90</td><td>3.50</td><td>0.01</td><td>0.12</td><td>0.91</td><td>-11.11</td><td>0.10</td><td>0.67</td><td>0.75</td><td>12.99</td><td>0.61</td><td>0.59</td><td>0.90</td><td>3.50</td><td>0.52</td></tr><tr><td>I'</td><td>Turia</td><td>0.35</td><td>0.87</td><td>-25.62</td><td>0.71</td><td>0.27</td><td>0.75</td><td>-27.75</td><td>0.43</td><td>0.62</td><td>0.81</td><td>13.74</td><td>0.54</td><td>0.61</td><td>0.72</td><td>-10.63</td><td>0.41</td></tr></tbody></table>

### 3.1.2. Streamflow Evaluation of Additional Calibration (AC)

The monthly streamflow simulated with the Additional Calibration AC (see Figure 3 for more information about the procedure) are presented in Figure 4 and corresponding statistics in Table 3. The sensitivity parameters are presented in Supplementary Materials Table S2 (parameters sensitivity ranking) and Table S3 (sensitivity analysis). The most sensitive parameters differ depending on the SUDOE regions. However, the base flow alpha factor, the groundwater delay and the soil evaporation compensation factor are the most sensitive over the entire territory. According to Moriasi et al. [13,74] and Gupta et al. [76] standards, the statistical performance over the entire territory is satisfactory to good for calibration and validation periods with respectively, with a range of NSE of 0.52 (0.22; 0.72) and 0.46 (-0.01; 0.78), a R² of 0.78 (0.53; 0.94) and 0.76 (0.34; 0.94), a PBIAS of ±12.52 (-38.91; 13.74) and ±11.25 (-28.21; 6.00), and finally a KGE of 0.53 (0.15; 0.73) and 0.52 (0.27; 0.82).

Over the SUDOE main watersheds (Figure 2), AC simulated and observed annual irrigation volumes are rather similar (16,917 Hm³ and 17,537 Hm³, respectively). The regional comparison of the simulated volumes of irrigation from the AC calibration in each Spanish region (Figure 5) indicates also a good regional distribution of irrigated water in the model with R² (0.87) and PBIAS (10%).

The calibrated irrigation volumes were obtained by offsetting volumes per hectare, as the model overestimates the irrigated area because of the initial assumption (each type of land use was considered as totally irrigated or non-irrigated—Section 2.4.3). For this reason, the absolute irrigated volume spread by fields by hectares is partially biased with less water applied and a loss of biomass production. However, in order to compensate for the lack of biomass, the plant biomass production parameters were modified (cf. Supplementary Materials Table S1).
<figure><img src="image_8.png" /><figcaption>**Figure 5.** Simulated additional calibration (AC) versus observed annual water volume (Hm³) used in irrigation for each Spanish region (inter-annual averages from 2000 to 2010).</figcaption></figure>

### 3.1.3. Comparison of Streamflow Simulations between Conventional Calibration (CC) and Additional Calibration (AC)

CC and AC streamflow dynamics are similar (Figure 4) whereas the calibration procedures were different—AC procedure is enhanced by calibrating anthropogenic activities such as irrigation volume and crop yield. The parameters used to calibrate irrigation volume and crop yield for AC procedure are presented in Supplementary Materials Table S1. The sensitive parameters ranking according to the significance of the sensitivity (p-value) is introduced on the Supplementary Materials in Table S2 with the most sensitive parameters ranked first. CC and AC calibration are different sensitive parameters for streamflow calibration. However, similar sensitive parameters are found as the base flow alpha factor or the groundwater delay. Previous studies [11,16,92] show similar sensitivity for both of these parameters indicating a direct relationship between the surface runoff and other compartments parameters.

Water balances between CC and AC calibrations are different. The difference can be important ranging from -114% to +39% between compartments for surface runoff (-114%; 28%), lateral flow (-29%; 31%), water yield (-27%; 39%), evapotranspiration (-51%; 4%), and sediment yield (-100%; 10%). The strongest differences are presented during low flow periods with a difference of -18% on average whereas the lowest differences occur during high flow (-9% on average). It also seems that there are fewer differences for watersheds with low anthropogenic pressures (difference of -11%) whereas the difference is more important for watersheds highly anthropized (difference of -39%). However, no data from these hydrological components (lateral flow, evapotranspiration, sediment yield) with a sufficient temporal and spatial resolution can allow the validation of one or the other calibration process regarding the distribution of water fluxes inside the water cycle. Additional validation was conducted with the comparisons of water irrigated and crop yields simulations and national datasets values.
The total irrigated volume of AC (20,400 Hm³) is close to observed values (21,100 Hm³) (R² = 0.87, p-value < 0.001) compared to the CC simulations with the observed values (R² = 0.37, p-value < 0.001) which overestimates irrigated volumes by 36% (27,800 Hm³). Irrigation volume has a substantial impact on water repartition over the watershed and, consequently, on the water balance between compartments (surface runoff (-114%; 28%) and evapotranspiration (-51%; 4%)). In AC calibration, irrigation volume is taken from dams and rivers whereas, in CC calibration, water use for irrigation comes from outside the watershed. In CC water is not limiting and water management is not constrained. The choices of irrigation parameterization have an impact on water balance and hydrological response with the simulated watershed. This study highlights the importance for the SWAT model users to carefully verify the specification and boundaries of crop management and irrigation routines, in order to avoid an over or underestimation of the water volume. A suggestion is to compare the simulated crop management techniques with national statistic values to evaluate the truthfulness of the simulation.

### 3.1.4. Water Yield Validation with Observations and Model Outputs

A comparison of the water yield obtained from the different calibration methodologies has been conducted in order to validate the hydrological simulation on different sources of data and not just streamflow records. The water yield obtained from AC simulations and CC simulations was first compared with observed annual water volume (Figure 6A) in order to compare the performance of each calibration methodology. In a second step, the CC and AC projects once calibrated were validated by comparing water volumes under "natural" conditions with water volumes obtained with SIMPA under "natural" conditions as well (Figure 6B).

<figure><img src="image_9.png" /><figcaption>**Figure 6.** Simulated CC (conventional calibration) and AC (additional calibration) versus (A) observed annual water volume and versus (B) SIMPA simulation for the main watersheds of the SUDOÉ territory (inter-annual averages from 2000 to 2010).</figcaption></figure>

Figure 6A shows satisfactory regressions between the annual water volumes observed and simulated with CC ($R^2$ = 0.91, slope = 1.13) and AC calibrations ($R^2$ = 0.94, slope = 1.07). The good simulation of annual water yield supports the analysis of an acceptable simulation of discharge as described above. High $R^2$ was found between SWAT simulated "natural" interannual water volume (the simulation where all human pressures have been removed) and SIMPA simulated water volume (0.88 and 0.91, respectively for CC and AC; Figure 6B). The regression is good for both models but the slopes are different: 0.74 for CC and 0.99 for AC. The slope difference could be explained by the difference in irrigation calibration. To achieve a good calibration, the CC calibration was set
up to provide irrigation with water originating from outside the watershed inducing an increase of surface runoff that is compensated with the increase of infiltration into groundwater. For the AC regressions ($R^2$ > 0.9), linear slopes are close to 1 indicating that water volume in SWAT is not under- or overestimated. Water volumes seem correctly distributed between the main watersheds of the SUDOE territory for both “natural” and “anthropogenic” versions of the AC calibration. The double comparison of water volume (SIMPA and observation) reinforces the validation of AC simulations.

Monthly water volumes were then compared by subsystems (Figure 1B) in order to validate the spatial and temporal variations of the water volume on a smaller spatial scale. The $R^2$ between SIMPA and AC “natural” simulations is ranging between 0.84 and 0.98 (Figure 7), ascertaining a proper distribution of water volume within each main watershed.

<figure><img src="image_10.png" /><figcaption>**Figure 7.** “Natural” water volume for subsystem simulated with SWAT (AC; additional calibration) versus “natural” water volume simulated with SIMPA (Figure 1B) of the main Spanish watersheds (inter-annual averages from 2000 to 2010).</figcaption></figure>

Figure 8 shows the interannual monthly dynamics of “natural” water volume for the main watersheds between SIMPA, CC, and AC simulations. Two groups of watersheds can be defined: One group with high similarity in term of temporal dynamics of monthly water volumes (Ebro, Duero, Guadiana, Guadalquivir, Tajo) and another group with a significant difference between SWAT and SIMPA water volumes from June to December (Jucar and Segura). These two groups might be related to the level of anthropogenic pressures on the watershed. The Jucar and Segura watersheds, the most important ones in the South East of Spain, characterized by a warm and dry climate [93] with low flow support based on water imported from other watersheds, thus highly anthropized. The SIMPA model (the reference in this study), which modeled upstream streamflow with low pressures, seems to support water volumes almost stable throughout the year whereas SWAT simulation shows a severe decrease of water during low flow (Figure 8). Overall, most of the watersheds have the same trend between models hence decreasing the water balance uncertainties. Model comparison, with models that had attained a level of familiarity and trust by stakeholders, is efficient when dataset availability is limited and show promise to evaluate scenarios and decision-making [94,95].
<figure><img src="image_11.png" /><figcaption>**Figure 8.** Inter-annual monthly dynamics of “natural” water volumes in the main SUDOE watersheds simulated with SIMPA (red line), SWAT additional calibration (SWAT AC, blue line) and SWAT conventional calibration (SWAT CC, dotted dark line). The color ribbons (blue for AC and gray for CC) show the 95% prediction uncertainty intervals resulting from the best SWAT-CUP simulation. The letters (A,C,E–K) indicate the watershed localization over SUDOE territory (Figure 2) (inter-annual averages from 2000 to 2010). There is no SIMPA simulation for Garonne and Aude watersheds.</figcaption></figure>

## 3.2. Water Quality Evaluation

To evaluate water quality, sediment and nitrate loads are commonly used as indicators of the particulate and dissolved in-stream elements, respectively. Modeled sediment and nitrate load time series are compared to observations as well as the impact of cultural practices on CC and AC simulations.
### 3.2.1. Sediment Loads

The monthly sediment loads from LOADEST and both SWAT versions are shown in Figure 9. The mean $R^2$ of monthly sediment load is 0.72 ranging from 0.17 to 0.95. Sediment load is difficult to calibrate especially during flood events with greater deviations but the simulation gives a good understanding of the global trend [20,90].

<figure><img src="image_12.png" /><figcaption>**Figure 9.** Monthly variation of sediment loads from 2000 to 2010 (calibration period: 2000–2005, validation period: 2006–2010) obtained with LOADEST and SWAT calibrations in the main SUDOE watersheds (AC: additional calibration; CC: conventional calibration). The letters (A,C,D,F,I–K) indicate the watershed localization over SUDOE territory (Figure 2).</figcaption></figure>

The sediment loads between CC and AC calibrations are similar for most of the watersheds such as Aude, Minho, Guadiana, and Segura (Figure 9). However, CC and AC calibrations do not have the same response in some other watersheds such as the outlet of Garonne, Ebro, and Jucar (Figure 9A,F,I) where the model tends to overestimate the sediment loads in AC calibration whereas these loads are underestimated with CC calibration. The Garonne, Ebro, and Jucar are all main agricultural watersheds. The sediment loads difference is due to the increase in water volume used in irrigation that increases the soil runoff and erosion coefficient. In the Garonne, Ebro, and Jucar watersheds, the main difference between CC and AC is the localization of the irrigation source. By default, the value is “outside the watershed” when the difference was changed. This irrigation source information could create a substantial difference in the streamflow in a watershed with low water levels and high irrigation such as Jucar watershed. The pedoclimatic conditions have indeed a substantial influence on
some parameters impacting sediment load. In fact, the weather, the lithology and the land use have an impact on soil erosion and some parameters such as USLE C factor or the management practice will have an impact on this erosion [47].

### 3.2.2. Nitrate Loads

The interpolation of nitrate concentration with LOADEST gives satisfactory results (Supplementary Materials Figures S2 and S3) according to Kim et al. [96] standards. LOADEST tends to overestimate extreme values as largely [97], but this work is focusing on the impact of global change on water resources and more precisely on the SUDOE water scarcity issues which mostly occur during low flow periods [33]. Other studies already compared the nitrate loads simulated by LOADEST and SWAT models [97] and showed that the LOADEST software is useful when you have a discrete and large range of values under different hydrological conditions.

The comparison of nitrate loads between SWAT versions and LOADEST at the main river outlets shows an overall good agreement with an $R^2$ of 0.67 (0.55; 0.82) and an NSE of 0.51 (−0.36; 0.67) for the calibration period and with an $R^2$ of 0.65 (0.34; 0.93) and an NSE of 0.42 (−0.36; 0.91) for the validation period (Figure 10). The differences between AC and CC calibrations are not significant (*p*-value > 0.29) but with respect to the absolute performance value, nitrate loads of the CC calibration have lower statistics than the AC calibration with a mean $R^2$ of 0.50 and 0.70, respectively. For Jucar and Guadiana watersheds, nitrate loads of CC calibration (Figure 10) are lower than AC which could be explained by a decrease of nitrate runoff correlated to an increase of biomass production produced by the non-limiting irrigation. There are two actions that influence this response. The first is auto-irrigation depending on the plant demand (the production is not limited) and the second is fertilization. Biomass production is high thanks to the presence of fertilizer, applied in the most effective way, decreasing the soluble soil nitrate concentration which leads to a decrease in the nitrate runoff. For other watersheds such as Garonne, Aude, or Segura, nitrate runoff is more important during CC calibration because irrigation is triggered according to soil water saturation, causing excessive irrigation in these basins, which leads to nitrate runoff after fertilizer application. The comparison between observed and simulated crop yields is thus essential to determine nitrate sources and to validate the model. Jégo et al. [98] indicate the importance of crop management and particularly for N application periods on nitrate runoff in the Garonne watershed. Hence, it is therefore important to carefully take into account crop management, as well as crop yield and fertilization in the model.

### 3.2.3. Yield and Fertilization

Observed and AC simulated interannual crop yields over the basin show a good adjustment according to the average PBIAS indicator of +8% (−78; +40) (Figure 11). The average crop yield over the territory was 7.4 $t \ ha^{-1}$ for the AC calibration but only 2.3 $t \ ha^{-1}$ for CC. Two hypotheses might explain the CC crop yield underestimations: (1) The fertilization is too low, (2) the irrigation is only triggered when the plant is under stress which can slow down plant growth.
<figure><img src="image_13.png" /><figcaption>**Figure 10.** Average variation of daily nitrate loads from 2000 to 2010 (calibration period: 2000–2005: validation period: 2006–2010 simulated at monthly scale) obtained with Load Estimator model (LOADEST) and SWAT versions in the main SUDOE watersheds (AC: additional calibration; CC: conventional calibration). The letters (A,C,D,F,G,I–K) indicate the watershed localization over SUDOE territory (Figure 2).</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>**Figure 11.** Boxplot of annual national and AC simulated crop yield values (AC: additional calibration (A) by country; (B) by crops. AC simulated crop yields are the interannual average from 2000 to 2010.)</figcaption></figure>
### 3.3. Variability of the Model Through the Climate Change Scenario

Table 4 represents the percentage of change of each variable for both calibrations with a –20% decrease in precipitation. A negative impact on variables is expected with a diminishing rainfall such as a reduction in the irrigated volume and crop yields [99]. CC and AC water yields have the same trend with a decrease close to –40%. A streamflow decrease is also highlighted by Stahl et al. [100]. However, the crop management variables are different between CC and AC because CC crop management is not calibrated. Under drier conditions, plant growth will be coerced by the lack of water in the watershed in the AC calibration (–35%) whereas irrigated water is increasing in the CC calibration (+5%). The same opposite trends occurred with crop yield simulation with respectively +23% and –26% for CC and AC calibration. The irrigation source used in the CC is located outside of the watershed. Under such condition, the auto-irrigation is not limited by the water volume available in the river and the plant growth are not impacted by the drought period. That is the reason why the model simulated a light increase of plant growth with a decrease of rainfall. AC outputs are the most realistic calibration considering that plant growth, crop yield and irrigation volume will be coerced by the lack of water under a drier climate. This study strongly recommends calibrating crop management in the water scarcity context in order to have a proper allocation of water over the watershed, especially for studies assessing the potential impacts of climate change or human practices changes. The question is: Can we have an indicator to help us in the choice of calibration/validation procedure for one watershed?

<table><caption>**Table 4.** Percentage of variable differences for conventional calibration (CC) and additional calibration (AC) between the regular precipitation and a decrease of 20% of rainfall from 2000 to 2010. A positive percentage represents an increase of the variable whereas a negative value means a decrease.</caption><thead><tr><td>Variables</td><td>Differences of CC</td><td>Differences of AC</td></tr></thead><tbody><tr><td>Precipitation</td><td>-20%</td><td>-20%</td></tr><tr><td>Water yield</td><td>-44%</td><td>-41%</td></tr><tr><td>Fertilization</td><td>17%</td><td>-9%</td></tr><tr><td>Biomass production</td><td>23%</td><td>-26%</td></tr><tr><td>Crop yield average</td><td>33%</td><td>-31%</td></tr><tr><td>Irrigation volume</td><td>5%</td><td>-35%</td></tr></tbody></table>

### 3.4. Index of Hydrological Alteration (IHA) Helping in the Choice of Calibration/Validation Procedure between CC and AC

In order to decide which calibration/validation procedure is the most appropriate, this study analyzed the variations of the IHA according to the deviation of KGE (ΔKGE). Figure 12 shows a negative correlation between the IHA and ΔKGE for the 47 gauging stations representing the entire range of subbasins over the SUDOE territory. Two groups can be distinguished: (1) Subbasins with IHA above 80% and (2) subbasins under 80%. With an IHA above 80%, the ΔKGE is lower than 0.1, the CC and AC calibrations have similar trends so it might not be useful to calibrate the human activities as dams and crop management. When the IHA is under 80%, the ΔKGE variations are high (0; 0.65)—meaning that in some of these cases, the addition and the calibration of human activities, such as dams and crop management, is necessary to obtain a representative simulation of the reality. Our recommendation will be to apply the AC procedure in the case of when the IHA is under 80%. The IHA index can become a tool to determine the calibration/validation procedure. It is available in the south-west region but could be transferable to regions where water resources are impacted by agriculture and water management (dams, channels).
<figure><img src="image_15.png" /><figcaption>**Figure 12.** ΔKGE between conventional (CC) and additional calibration (AC) for 47 gauging stations face to IHA (Hydrological Alteration Index) from 2000 to 2010.</figcaption></figure>

This study explores an alternative methodology to model water quantity and quality resources at multiple scales, for the first time, by coupling various data sources (agricultural, urban effluents, quality and quantity data) at the SUDOE scale—that has been never approached [101–103]. Previous studies made analysis on uncalibrated simulations [9] or calibrated only on streamflow [11]. For hydrological models, only the streamflow and the water balance are usually analyzed based on gauging station records and previous studies [8]. For agronomic models, biomass production and crop yield help to calibrate the model [104]. This study combines hydrological and agronomic calibration approaches in order to improve the model estimation thanks to multiple datasets with different scales and time resolutions. In this study, a suggestion is made: Calibrating the model with different data sources, easily available coming from national open-access sources such as national statistics department or water agency reports and samplings. The multiple and diversified databases used showed that the water balance and human activities were better represented in AC calibration. The multi-watershed model is complex. Sometimes a trade-off of one parameter needs to be found because watershed responses to the parameter could be opposite. The need to cross data for calibrating and validating models is underlined. This study should not be reduced to the SWAT model only but it could be extrapolated to every hydro-agronomic model. Outputs from calibrated models come along with model uncertainties and assumptions. It is up to the analysts to quantify model prediction uncertainty and calibrate with their own field expertise of the watershed. In order to establish this methodology, running a soft calibration first and then launching an optimization procedure is advised. The soft calibration is a manual calibration based on watershed knowledge (characteristics, specificities). The optimization procedure is allowed to refine the parameters calibrated during the soft calibration, as a statistical calibration, for example, SWAT-CUP for SWAT [69,71].

This paper provides some tools and feedbacks on the SWAT model set up in order (1) to understand the importance of datasets choices in the calibration and validation steps, (2) to help future users in the large scale set up and (3) to choose the calibration/validation procedure. So far, the AC calibration is thus the most appropriate simulation to run further research on ecosystem services evaluations, climate change impact or direct human activities impact (e.g., the effect of implementing new water policy).

### 3.5. Limitations and Uncertainties of This Study

We calculated the percentage of rivers of Table 3 the statistics criteria R², NSE, PBIAS and KGE. According to Moriasi [13,74], we considered the performances of the model succeed when R², is higher to 0.6, NSE to 0.5, and PBIAS ranged between >15% or <-15%. According to the Knoben criteria [77],
we considered the performance succeed when the KGE is higher to -0.41. Based on these criteria, the hydrology performances of the model for both calibration and validation period and both conventional and additional calibration are succeeded 100% for KGE, 90% for R², 67% for PBIAS, and 31% for NSE. The performances obtained in this work for some rivers highly anthropized in the South of Europe are similar or better than the ones published by Abbaspour et al. [9]. Previous studies tried to compare observed and modeled continental river discharge and nitrate loads [9,105]. These studies point out the difficulties to calibrate and validate hydrological processes at the continental scale and are rather more complex for nitrate loads. In fact, our study model had to simplify the model set up and the input because of the lack of data related to anthropization. The quality of the outputs is significantly impacted by the degree of simplification, the input data availability and uncertainty [106], the temporal and spatial degree of discretization [107] which also impacts the quality of calibration and validation procedure. In the context of high anthropized systems, the performances of the calibration are related to water resources management (dams, canals, irrigation practices). All these data set should be at the scale of daily time step but unfortunately for dams, canals, and crop practices, we have data set at an annual time scale. In this context, the input data increase the uncertainties of the output. These uncertainties also impact the nitrate and sediment loads and the associated processes. The SWAT model used simplify equations [43,47] for quantifying the sediment and nitrate loads into the river. Therefore, these simplifications can cause deviations for sediment routing and nitrogen modeling. In fact, the uncertainty in model inputs influenced the N soil input such as the fertilizer use or the allocation of animal manure (that we do not consider in this study). The uncertainties in the exchange fractions of different types of N are still not completely known. Improving N leaching and runoff fraction should receive more attention in future studies in order to improve the nitrogen cycle modeling.

In order to optimize the calibration procedure, in this paper, we proposed a method to evaluate the level of the anthropization by the IHA index for large systems and we show that the model performance is directly related to this index (low IHA, high pressure, and low-performance model). In this context, the application of the Moriasi statistical thresholds is not adapted. Finally, our study has some limitations and high uncertainties which should be taking into consideration for further studies.

## 4. Conclusions

This study shows for the first time the evaluation of a guideline for calibrating multi-watersheds hydrological models as a complexity according to a level of anthropization in the context of global changes. The IAH is a good index to determine the complexity level of calibration procedure. To determine the most appropriate calibration/validation procedures, the global hydrological alteration index, IHA, demonstrated its feasibility: IHA values under 80% indicate that a complex procedure of calibration/validation would yield better results. Many applications of the model should be done in further studies. An assessment of ecosystem services and ecological functions could be interesting as well as the economic evaluations of these services such as the economic cost of maintaining good water quality. The model might also help to conduct river and land use planning studies in order to help stakeholders. Finally, the model can be used to conduct climate and land-use change impact studies.

**Supplementary Materials:** The following are available online at http://www.mdpi.com/2073-4441/12/1/115/s1, Figure S1: Comparison between observed nitrogen loads (coming from the UWWTP databases) and calculated one with the Zessner and Lindtner (2005) method, Figure S2: Sediment load regression between observation and LOADEST simulation at daily time step in kg d⁻¹, Figure S3: Nitrate load regression between observation and LOADEST simulation at daily time step in kg d⁻¹, Table S1: Parameters used to calibrate water quantity and quality in the SUDOE territory. The parameters to calibrate irrigated volume and crop yield were used only for the AC calibration, Table S1: Parameters sensitivity ranking for streamflow model outputs for each part of SUDOE territory for CC and AC models, Table S2: Sensitivity analysis for CC streamflow models outputs, Table S3: Sensitivity analysis for AC streamflow models outputs.

**Author Contributions:** Conceptualization, R.C., M.R., S.S. and J.M.S.-P.; Data curation, R.C., J.P.-A., Y.G., L.R., E.N., M.S.-C., J.L.L.-C., J.J.G.P. and J.M.B.; Formal analysis, R.C., J.P.-A. and L.R.; Funding acquisition, S.S. and J.M.S.-P.; Investigation, R.C., M.R., M.M., J.L.L.-C. and J.M.S.-P.; Methodology, R.C., M.R. and S.S.; Project administration, S.S., J.P.-A., E.N. and J.M.S.-P.; Resources, J.M.S.-P.; Supervision, S.S. and J.M.S.-P.; Visualization,
R.C.; Writing—original draft, R.C.; Writing—review and editing, M.R., S.S., Y.G., M.M., E.N., J.L.L.-C., J.J.G.P. and J.M.S.-P. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Acknowledgments:** This project was supported by the Université Toulouse III Paul Sabatier and by the Ministry of Higher Education and Research. This study was performed as part of the EU Interreg SUDOE IVB program (AGUAMOD—SOE1/P5/F0026 project, http://www.aguamod-sudoe.eu) and also supported by ERDF. We want to thank all the partners of this program that shared data and helped with their expertise on irrigation, dam and crop yield management and streamflow calibration: l'Agence de l'eau Adour-Garonne, la Compagnie d'Aménagement des Coteaux de Gascogne, el Universidad Politecnica de Valencia, la Banque Hydro, l'Electricité De France (EDF) (the project REGARD-RTRA/STAE), le Syndicat Mixte d'Etude et d'Aménagement de la Garonne, l'Institute National de la Statistique et des Etudes Economiques, el Instituto Nacional de Estatística de Portugal y de España, el Consejo Superior de Investigaciones Científicas, El Ministerio para la Transición Ecológica and el Ministerio de Agricultura, Pesca y Alimentación of Spain.

**Conflicts of Interest:** The authors declare no conflict of interest.

## References

1. Eliasson, J. The rising pressure of global water shortages. *Nat. News* **2015**, *517*, 6. [CrossRef]

2. Mekonnen, M.M.; Hoekstra, A.Y. Four billion people facing severe water scarcity. *Sci. Adv.* **2016**, *2*, e1500323. [CrossRef]

3. Vörösmarty, C.J.; McIntyre, P.B.; Gessner, M.O.; Dudgeon, D.; Prusevich, A.; Green, P.; Glidden, S.; Bunn, S.E.; Sullivan, C.A.; Liermann, C.R.; et al. Global threats to human water security and river biodiversity. *Nature* **2010**, *467*, 555–561. [CrossRef] [PubMed]

4. European Environment Agency (EEA). EEA Report Water Use and Environmental Pressures; EEA: Copenhagen, Denmark, 2018.

5. Raimonet, M.; Thieu, V.; Silvestre, M.; Oudin, L.; Rabouille, C.; Vautard, R.; Garnier, J. Landward Perspective of Coastal Eutrophication Potential Under Future Climate Change: The Seine River Case (France). *Front. Mar. Sci.* **2018**. [CrossRef]

6. Roudier, P.; Andersson, J.C.M.; Donnelly, C.; Feyen, L.; Greuell, W.; Ludwig, F. Projections of future floods and hydrological droughts in Europe under a +2 °C global warming. *Clim. Chang.* **2016**, *135*, 341–355. [CrossRef]

7. Friesen, J.; Rodríguez-Sinobas, L.; Barceló, D. *Advanced Tools for Integrated Water Resources Management*; Academic Press: Cambridge, MA, USA, 2018; ISBN 978-0-12-814300-1.

8. Fu, B.; Merritt, W.S.; Croke, B.F.W.; Weber, T.R.; Jakeman, A.J. A review of catchment-scale water quality and erosion models and a synthesis of future prospects. *Environ. Model. Softw.* **2019**, *114*, 75–97. [CrossRef]

9. Abbaspour, K.C.; Rouholahnejad, E.; Vaghefi, S.; Srinivasan, R.; Yang, H.; Kløve, B. A continental-scale hydrology and water quality model for Europe: Calibration and uncertainty of a high-resolution large-scale SWAT model. *J. Hydrol.* **2015**, *524*, 733–752. [CrossRef]

10. Falkenmark, M.; Rockstrom, J.; Rockström, J. *Balancing Water for Humans and Nature: The New Approach in Ecosystems*; Earthscan: London, UK, 2004; ISBN 978-1-85383-926-9.

11. Grusson, Y.; Anctil, F.; Sauvage, S.; Sánchez Pérez, J.M. Coevolution of Hydrological Cycle Components under Climate Change: The Case of the Garonne River in France. *Water* **2018**, *10*, 1870. [CrossRef]

12. Biancamaria, S.; Mballo, M.; Le Moigne, P.; Sánchez Pérez, J.M.; Espitalier-Noël, G.; Grusson, Y.; Cakir, R.; Häfliger, V.; Barathieu, F.; Trasmonte, M.; et al. Total water storage variability from GRACE mission and hydrological models for a 50,000 km² temperate watershed: The Garonne River basin (France). *J. Hydrol. Reg. Stud.* **2019**, *24*, 100609. [CrossRef]

13. Moriasi, D.; Gitau, M.W.; Pai, N.; Daggupati, P. Hydrologic and Water Quality Models: Performance Measures and Evaluation Criteria. *Trans. ASABE* **2015**, *58*, 1763–1785.

14. Holzworth, D.P.; Huth, N.I.; de Voil, P.G.; Zurcher, E.J.; Herrmann, N.I.; McLean, G.; Chenu, K.; van Oosterom, E.J.; Snow, V.; Murphy, C.; et al. APSIM—Evolution towards a new generation of agricultural systems simulation. *Environ. Model. Softw.* **2014**, *62*, 327–350. [CrossRef]

15. Wade, A.J.; Durand, P.; Beaujouan, V.; Wessel, W.W.; Raat, K.J.; Whitehead, P.G.; Butterfield, D.; Rankinen, K.; Lepisto, A. A nitrogen model for European catchments: INCA, new model structure and equations. *Hydrol. Earth Syst. Sci.* **2002**, *6*, 559–582. [CrossRef]
16. Arnold, J.G.; Srinivasan, R.; Muttiah, R.S.; Williams, J.R. Large Area Hydrologic Modeling and Assessment Part I: Model Development. *JAWRA J. Am. Water Resour. Assoc.* **1998**, *34*, 73–89. [CrossRef]

17. Ruelland, D.; Billen, G.; Brunstein, D.; Garnier, J. SENEQUE: A multi-scaling GIS interface to the Riverstrahler model of the biogeochemical functioning of river systems. *Sci. Total Environ.* **2007**, *375*, 257–273. [CrossRef] [PubMed]

18. Sparrow, A.D.; Friedel, M.H.; Smith, D.M.S. A landscape-scale model of shrub and herbage dynamics in Central Australia, validated by satellite data. *Ecol. Model.* **1997**, *97*, 197–216. [CrossRef]

19. Valigura, R.A.; Alexander, R.B.; Castro, M.S.; Meyers, T.P.; Paerl, H.W.; Turner, R.E.; Stacey, P.E. *Nitrogen Loading in Coastal Water Bodies: An Atmospheric Perspective*; American Geophysical Union: Washington, DC, USA, 2001; ISBN 978-0-87590-271-5.

20. Oeurng, C.; Sauvage, S.; Sánchez-Pérez, J.-M. Assessment of hydrology, sediment and particulate organic carbon yield in a large agricultural catchment using the SWAT model. *J. Hydrol.* **2011**, *401*, 145–153. [CrossRef]

21. Volk, M.; Sauvage, S.; Sánchez-Pérez, J.-M.; Biegerd, K. Development and applications of the SWAT model to support sustainable river basin management on different scales. *Sustain. Water Qual. Ecol.* **2016**, *8*, 1–3. [CrossRef]

22. Gassman, P.W.; Sadeghi, A.M.; Srinivasan, R. Applications of the SWAT Model Special Section: Overview and Insights. *J. Environ. Qual.* **2014**, *43*, 1–8. [CrossRef]

23. Brown, L.C.; Barnwell, T.O. *The Enhanced Stream Water Quality Models QUAL2E and QUAL2E-UNCAS: Documentation and User Model*; Environmental Research Laboratory, Office of Research and Development, U.S. Environmental Protection Agency: Washington, DC, USA, 1987.

24. Sharpley, A.N.; Williams, J.R. EPIC-Erosion/Productivity Impact Calculator: 1, Model Documentation. *USDA Agric. Res. Serv. Tech. Bull.* **1990**, *235*, 1759.

25. Arnold, J.G.; Srinivasan, R.; Muttiah, R.S.; Allen, P.M. Continental scale simulation of the hydrologic balance. *JAWRA J. Am. Water Resour. Assoc.* **1999**, *35*, 1037–1051. [CrossRef]

26. Kinerson, R.S.; Kittle, J.L.; Duda, P.B. BASINS: Better Assessment Science Integrating Point and Nonpoint Sources. In *Decision Support Systems for Risk-Based Management of Contaminated Sites*; Marcomini, A., Suter II, G.W., Critto, A., Eds.; Springer US: Boston, MA, USA, 2009; pp. 1–24. ISBN 978-0-387-09722-0.

27. Zalasiewicz, J.; Waters, C.; Summerhayes, C.; Williams, M. The Anthropocene. *Geol. Today* **2018**, *34*, 177–181. [CrossRef]

28. Lévêque, C. *La Mémoire des Fleuves et des Rivières—Lévêque*; Ulmer: Paris, France, 2019; ISBN 978-2-84138-968-1.

29. Volk, M.; Bosch, D.; Nangia, V.; Narasimhan, B. SWAT: Agricultural water and nonpoint source pollution management at a watershed scale. *Agric. Water Manag.* **2016**, *175*, 1–3. [CrossRef]

30. Volk, M. Modelling ecosystem services: Current approaches, challenges and perspectives. *Sustain. Water Qual. Ecol.* **2015**, *5*, 1–2. [CrossRef]

31. Grizzetti, B.; Pistocchi, A.; Liquete, C.; Udias, A.; Bouraoui, F.; van de Bund, W. Human pressures and ecological status of European rivers. *Sci. Rep.* **2017**, *7*, 205. [CrossRef]

32. Richter, B.D.; Baumgartner, J.V.; Powell, J.; Braun, D.P. A Method for Assessing Hydrologic Alteration within Ecosystems. *Conserv. Biol.* **1996**, *10*, 1163–1174. [CrossRef]

33. EEA Report European Waters—Assessment of Status and Pressures 2018. Available online: https://www.eea.europa.eu/publications/state-of-water (accessed on 22 March 2019).

34. Eurostat. The Statistical Office of the European Union Water use by NUTS 2 regions—Eurostat. Available online: https://ec.europa.eu/eurostat/web/products-datasets/-/env Watuse r2 (accessed on 14 May 2019).

35. Romero, E.; Garnier, J.; Lassaletta, L.; Billen, G.; Le Gendre, R.; Riou, P.; Cugier, P. Large-scale patterns of river inputs in southwestern Europe: Seasonal and interannual variations and potential eutrophication effects at the coastal zone. *Biogeochemistry* **2013**, *113*, 481–505. [CrossRef]

36. Gao, Y.; Vogel, R.M.; Kroll, C.N.; Poff, N.L.; Olden, J.D. Development of representative indicators of hydrologic alteration. *J. Hydrol.* **2009**, *374*, 136–147. [CrossRef]

37. Richter, B.D.; Mathews, R.; Harrison, D.L.; Wigington, R. Ecologically Sustainable Water Management: Managing River Flows for Ecological Integrity. *Ecol. Appl.* **2003**, *13*, 206–224. [CrossRef]

38. Richter, B.D.; Baumgartner, J.V.; Braun, D.P.; Powell, J. A spatial assessment of hydrologic alteration within a river network. *Regul. Rivers Res. Manag.* **1998**, *14*, 329–340. [CrossRef]

39. Büttner, G.; Kosztra, B. *Addendum to CLC2006 Technical Guidelines*; Version 2; EEA: Copenhagen, Denmark, 2014.
40. Rivas-Martínez, S.; Penas, A.; Díaz, T.E. Bioclimatic map of Europe, bioclimates, University of León, E-24071, Spain. 2004. Available online: https://www.globalbioclimatics.org/form/maps.htm (accessed on 29 December 2019).

41. AQUASTAT. Available online: http://www.fao.org/nr/water/aquastat/data/query/results.html (accessed on 3 September 2019).

42. USDA Soil Conservation Service Hydrology Section 4. In *National Engineering Handbook*; Natural Resources Conservation Service: Washington, DC, USA, 1972; Chapter 4–10.

43. Cunge, J.A. On the Subject of a Flood Propagation Computation Method (Musklingum Method). *J. Hydraul. Res.* **1969**, *7*, 205–230. [CrossRef]

44. Hargreaves, G.H.; Zohrab, A. Samani Reference Crop Evapotranspiration from Temperature. *Appl. Eng. Agric.* **1985**, *1*, 96–99. [CrossRef]

45. Hargreaves, G.H.; Allen Richard, G. History and Evaluation of Hargreaves Evapotranspiration Equation. *J. Irrig. Drain. Eng.* **2003**, *129*, 53–63. [CrossRef]

46. Williams, J.R. Sediment-yield prediction with Universal Equation using runoff energy factor. *Present Prospect. Technol. Predict. Sediment Yield Sources* **1975**, *40*, 244–252.

47. Neitsch, S.L.; Arnold, J.G.; Kiniry, J.R.; Williams, J.R. *Soil and Water Assessment Tool: Theoretical Documentation Version 2009*; Technical Report No. 406; Texas Water Resources Institute: Temple, TX, USA, 2011.

48. Gassman, P.; Reyes, M.; Green, C.; Arnold, J. The Soil and Water Assessment Tool: Historical Development, Applications, and Future Research Directions. *Trans. ASABE* **2007**, *50*, 1211–1250. [CrossRef]

49. Battude, M.; Al Bitar, A.; Brut, A.; Tallec, T.; Huc, M.; Cros, J.; Weber, J.-J.; Lhuissier, L.; Simonneaux, V.; Demarez, V. Modeling water needs and total irrigation depths of maize crop in the south west of France using high spatial and temporal resolution satellite imagery. *Agric. Water Manag.* **2017**, *189*, 123–136. [CrossRef]

50. Gilmore, E.C.; Rogers, J.S. Heat Units as a Method of Measuring Maturity in Corn 1. *Agron. J.* **1958**, *50*, 611–615. [CrossRef]

51. Monteith, J.L.; Moss, C.J.; Cooke, G.W.; Pirie, N.W.; Bell, G.D.H. Climate and the efficiency of crop production in Britain. *Philos. Trans. R. Soc. Lond. B Biol. Sci.* **1977**, *281*, 277–294. [CrossRef]

52. Gesch, D.; Farr, T.; Slater, J.; Muller, J.-P.; Cook, S. New products from the shuttle radar topography mission. *EOS Trans. Am. Geophys. Union* **2006**, *87*, 174. [CrossRef]

53. Nachtergaele, F.O.; van Velthuizen, H.; Verelst, L.; Batjes, N.H.; Dijkshoorn, J.A.; van Engelen, V.W.P.; Fischer, G.; Jones, A.; Montanarella, L.; Petri, M.; et al. *Harmonized World Soil Database (Version 1.0)*; FAO: Rome, Italy, 2008.

54. Wieder, W.R.; Boehnert, J.; Bonan, G.B.; Langseth, M. Regridded Harmonized World Soil Database v1.2. ORNL DAAC: Oak Ridge, TN, USA, 2014.

55. Håggmark, L.; Ivarsson, K.-I.; Gollvik, S.; Olofsson, P.-O. Mesan, an operational mesoscale analysis system. *Tellus A* **2000**, *52*, 2–20. [CrossRef]

56. Landelius, T.; Dahlgren, P.; Gollvik, S.; Jansson, A.; Olsson, E. A high-resolution regional reanalysis for Europe. Part 2: 2D analysis of surface temperature, precipitation and wind. *Q. J. R. Meteorol. Soc.* **2016**, *142*, 2132–2142.

57. Raimonet, M.; Oudin, L.; Thieu, V.; Silvestre, M.; Vautard, R.; Rabouille, C.; Le Moigne, P. Evaluation of Gridded Meteorological Datasets for Hydrological Modeling. *J. Hydrometeorol.* **2017**, *18*, 3027–3041. [CrossRef]

58. Quintana-Seguí, P.; Le Moigne, P.; Durand, Y.; Martin, E.; Habets, F.; Baillon, M.; Canellas, C.; Franchisteguy, L.; Morel, S. Analysis of Near-Surface Atmospheric Variables: Validation of the SAFRAN Analysis over France. *J. Appl. Meteorol. Climatol.* **2008**, *47*, 92–107. [CrossRef]

59. Quintana Seguí, P.; Peral, C.; Turco, M.; Llasat, M.C.; Martin, E. Meteorological analysis systems in North-East Spain: Validation of SAFRAN and SPAN. *J. Environ. Inf.* **2016**. [CrossRef]

60. EEA Report UWWTD Data Sources. Available online: https://www.eea.europa.eu/themes/water/european-waters/water-use-and-environmental-pressures/uwwtd/uwwtd-data-sources (accessed on 27 November 2018).

61. Zessner, M.; Lindtner, S. Estimations of municipal point source pollution in the context of river basin management. *Water Sci. Technol.* **2005**, *52*, 175–182. [CrossRef] [PubMed]

62. Campardon, M.; Loubier, S.; Morardet, S. *The Current Situation of Irrigation in France and Its Evolutions Since 2000—Rapport Final Convention BSE—IRSTEA 2011–2012*; CemOA: Antony, France, 2012; p. 58.
63. Gleyses, G.; Rieu, T. *L'irrigation en France: État des lieux 2000 et évolution*; Editions Quae: Versailles, France, 2004; ISBN 978-2-85362-642-2.

64. Loubier, S.; Campardon, M.; Morardet, S. L'irrigation diminue-t-elle en France? Premiers enseignements du recensement agricole de 2010. *Sci. Eaux Territ.* **2013**, *12–19*. [CrossRef]

65. Ministry of the Presidency. *Ley 10/2001, de 5 de julio, del Plan Hidrológico Nacional*; BOE-A-2001-13042; Ministry of the Presidency: Madrid, Spain, 2001.

66. Ministerio de Medio Ambiente. *Plan Hidrológico Nacional*; Secretaría de Estado de Aguas y Costas: Madrid, Spain, 2002.

67. Hengl, T.; de Jesus, J.M.; MacMillan, R.A.; Batjes, N.H.; Heuvelink, G.B.M.; Ribeiro, E.; Samuel-Rosa, A.; Kempen, B.; Leenaars, J.G.B.; Walsh, M.G.; et al. SoilGrids1k—Global Soil Information Based on Automated Mapping. *PLoS ONE* **2014**, *9*, e105992. [CrossRef] [PubMed]

68. Roosevelt, E. *Land Husbandry—Components and Strategy*; FAO: Rome, Italy, 1996.

69. Abbaspour, K.C.; Vaghefi, S.A.; Srinivasan, R. A Guideline for Successful Calibration and Uncertainty Analysis for Soil and Water Assessment: A Review of Papers from the 2016 International SWAT Conference. *Water* **2018**, *10*, 6. [CrossRef]

70. Van Griensven, A.; Meixner, T.; Grunwald, S.; Bishop, T.; Diluzio, M.; Srinivasan, R. A global sensitivity analysis tool for the parameters of multi-variable catchment models. *J. Hydrol.* **2006**, *324*, 10–23. [CrossRef]

71. Arnold, J.G.; Moriasi, D.N.; Gassman, P.W.; Abbaspour, K.C.; White, M.J.; Srinivasan, R.; Santhi, C.; Harmel, R.D.; van Griensven, A.; Van Liew, M.W.; et al. SWAT: Model Use, Calibration, and Validation. *Trans. ASABE* **2012**, *55*, 1491–1508.

72. Khalid, K.; Ali, M.F.; Rahman, N.F.A.; Mispan, M.R.; Haron, S.H.; Othman, Z.; Bachok, M.F. Sensitivity Analysis in Watershed Model Using SUFI-2 Algorithm. *Procedia Eng.* **2016**, *162*, 441–447. [CrossRef]

73. Yang, J.; Reichert, P.; Abbaspour, K.C.; Xia, J.; Yang, H.M. Comparing uncertainty analysis techniques for a SWAT application to the Chaohe Basin in China. *J. Hydrol.* **2008**, *358*, 1–23. [CrossRef]

74. Moriasi, D.N.; Arnold, J.G.; Van Liew, M.W.; Bingner, R.L.; Harmel, R.D.; Veith, T.L. Model evaluation guidelines for systematic quantification of accuracy in watershed simulations. *Trans. ASABE* **2007**. [CrossRef]

75. Gupta, H.V.; Kling, H.; Yilmaz, K.K.; Martinez, G.F. Decomposition of the mean squared error and NSE performance criteria: Implications for improving hydrological modelling. *J. Hydrol.* **2009**, *377*, 80–91. [CrossRef]

76. Kling, H.; Fuchs, M.; Paulin, M. Runoff conditions in the upper Danube basin under an ensemble of climate change scenarios. *J. Hydrol.* **2012**, *424*, 264–277. [CrossRef]

77. Knoben, W.J.M.; Freer, J.E.; Woods, R.A. Technical note: Inherent benchmark or not? Comparing Nash–Sutcliffe and Kling–Gupta efficiency scores. *Hydrol. Earth Syst. Sci. Discuss.* **2019**, *23*, 1–7. [CrossRef]

78. Runkel, R.L.; Crawford, C.G.; Cohn, T.A. Load estimator (LOADEST): A FORTRAN program for estimating constituent loads in streams and rivers. In *Techniques and Methods*; USGS: Washington, DC, USA, 2004. [CrossRef]

79. Boithias, L.; Srinivasan, R.; Sauvage, S.; Macary, F.; Sánchez-Pérez, J.M. Daily Nitrate Losses: Implication on Long-Term River Quality in an Intensive Agricultural Catchment of Southwestern France. *J. Environ. Qual.* **2014**, *43*, 46–54. [CrossRef]

80. Meaurio, M.; Zabaleta, A.; Uriarte, J.A.; Srinivasan, R.; Antigüedad, I. Evaluation of SWAT models performance to simulate streamflow spatial origin. The case of a small forested watershed. *J. Hydrol.* **2015**, *525*, 326–334. [CrossRef]

81. Mittelstet, A.R.; Storm, D.E.; White, M.J. Using SWAT to enhance watershed-based plans to meet numeric water quality standards. *Sustain. Water Qual. Ecol.* **2016**, *7*, 5–21. [CrossRef]

82. Spearman, C. The Proof and Measurement of Association between Two Things. *Am. J. Psychol.* **1904**, *15*, 72–101. [CrossRef]

83. Alvarez, J.; Sanchez, A.; Quintas, L. SIMPA, a GRASS based Tool for Hydrological Studies. *Int. J. Geoinformatics* **2005**, *1*, 14.

84. Taguas, E.V.; Gómez, J.A.; Denisi, P.; Mateos, L. Modelling the Rainfall-Runoff Relationships in a Large Olive Orchard Catchment in Southern Spain. *Water Resour. Manag.* **2015**, *29*, 2361–2375. [CrossRef]
85. Vicente, D.J.; Rodríguez-Sinobas, L.; Garrote, L.; Sánchez, R. Chapter Five—Applying Water Accounting Methods Through Statistical Data and Simulation Models: The Duero Transboundary Watershed. In *Advances in Chemical Pollution, Environmental Management and Protection*; Friesen, J., Rodríguez-Sinobas, L., Eds.; Advanced Tools for Integrated Water Resources Management; Elsevier: Amsterdam, The Netherlands, 2018; Volume 3, pp. 115–146.

86. Vautard, R.; Gobiet, A.; Sobolowski, S.; Kjellström, E.; Stegehuis, A.; Watkiss, P.; Mendlik, T.; Landgren, O.; Nikulin, G.; Teichmann, C.; et al. The European climate under a 2 °C global warming. *Environ. Res. Lett.* **2014**, *9*, 034006. [CrossRef]

87. Conaway, G.L.; Garen, D.C.; Wood-Wward, D.E.; Hoeft, C.C. Chapter 20—Watershed Yield. In *National Engineering Handbook Hydrology Chapters—NRCS*; NEH Part 630; United States Department of Agriculture: Washington, DC, USA, 2009; p. 29.

88. Santa-Maria, C.M.; Fernández Yuste, J.A. *Índices de Alteración Hidrológica en RíoS—Manual de Referencia Metodológica*; Dirección General el Agua (Ministerio de Medio Ambiente), Mediante Convenio con el CEDEX: Madrid, Spain, 2008.

89. Sun, N.-Z.; Sun, A. *Model Calibration and Parameter Estimation: For Environmental and Water Resource Systems*; Springer: Berlin, Germany, 2015; ISBN 978-1-4939-2323-6.

90. Wei, X.; Sauvage, S.; Le, T.P.Q.; Ouillon, S.; Orange, D.; Vinh, V.D.; Sanchez-Perez, J.-M. A Modeling Approach to Diagnose the Impacts of Global Changes on Discharge and Suspended Sediment Concentration within the Red River Basin. *Water* **2019**, *11*, 958. [CrossRef]

91. Jeuffroy, M.-H.; Casadebaig, P.; Debaeke, P.; Loyce, C.; Meynard, J.-M. Agronomic model uses to predict cultivar performance in various environments and cropping systems. A review. *Agron. Sustain. Dev.* **2014**, *34*, 121–137. [CrossRef]

92. Uribe, N.; Corzo, G.; Quintero, M.; van Griensven, A.; Solomatine, D. Impact of conservation tillage on nitrogen and phosphorus runoff losses in a potato crop system in Fuquene watershed, Colombia. *Agric. Water Manag.* **2018**, *209*, 62–72. [CrossRef]

93. Estrela, T.; Pérez-Martin, M.A.; Vargas, E. Impacts of climate change on water resources in Spain. *Hydrol. Sci. J.* **2012**, *57*, 1154–1167. [CrossRef]

94. Robles-Morua, A.; Halvorsen, K.E.; Mayer, A.S.; Vivoni, E.R. Exploring the application of participatory modeling approaches in the Sonora River Basin, Mexico. *Environ. Model. Softw.* **2014**, *52*, 273–282. [CrossRef]

95. Starkl, M.; Brunner, N.; López, E.; Martínez-Ruiz, J.L. A planning-oriented sustainability assessment framework for peri-urban water management in developing countries. *Water Res.* **2013**, *47*, 7175–7183. [CrossRef]

96. Kim, J.; Lim, K.J.; Park, Y.S. Evaluation of Regression Models of LOADEST and Eight-Parameter Model for Nitrogen Load Estimations. *Water Air. Soil Pollut.* **2018**, *229*, 179. [CrossRef]

97. Cerro, I.; Antigüedad, I.; Srinavasan, R.; Sauvage, S.; Volk, M.; Sanchez-Perez, J.M. Simulating land management options to reduce nitrate pollution in an agricultural watershed dominated by an alluvial aquifer. *J. Environ. Qual.* **2014**, *43*, 67–74. [CrossRef]

98. Jégo, G.; Martínez, M.; Antigüedad, I.; Launay, M.; Sanchez-Pérez, J.M.; Justes, E. Evaluation of the impact of various agricultural practices on nitrate leaching under the root zone of potato and sugar beet using the STICS soil-crop model. *Sci. Total Environ.* **2008**, *394*, 207–221. [CrossRef]

99. Olesen, J.E.; Carter, T.R.; Díaz-Ambrona, C.H.; Fronzek, S.; Heidmann, T.; Hickler, T.; Holt, T.; Minguez, M.I.; Morales, P.; Palutikof, J.P.; et al. Uncertainties in projected impacts of climate change on European agriculture and terrestrial ecosystems based on scenarios from regional climate models. *Clim. Chang.* **2007**, *81*, 123–143. [CrossRef]

100. Stahl, K.; Hisdal, H.; Hannaford, J.; Tallaksen, L.; Van Lanen, H.; Sauquet, E.; Demuth, S.; Fendekova, M.; Jordar, J. Streamflow trends in Europe: Evidence from a dataset of near-natural catchments. *Hydrol. Earth Syst. Sci.* **2010**, *14*, 2367–2382. [CrossRef]

101. Jeppesen, E.; Kronvang, B.; Olesen, J.E.; Audet, J.; Søndergaard, M.; Hoffmann, C.C.; Andersen, H.E.; Lauridsen, T.L.; Liboriussen, L.; Larsen, S.E.; et al. Climate change effects on nitrogen loading from cultivated catchments in Europe: Implications for nitrogen retention, ecological state of lakes and adaptation. *Hydrobiologia* **2011**, *663*, 1–21. [CrossRef]
102. Peraza-Castro, M.; Ruiz-Romera, E.; Meaurio, M.; Sauvage, S.; Sánchez-Pérez, J.M. Modelling the impact of climate and land cover change on hydrology and water quality in a forest watershed in the Basque Country (Northern Spain). *Ecol. Eng.* **2018**, *122*, 315-326. [CrossRef]

103. Pesce, M.; Critto, A.; Torresan, S.; Giubilato, E.; Santini, M.; Zirino, A.; Ouyang, W.; Marcomini, A. Modelling climate change impacts on nutrients and primary production in coastal waters. *Sci. Total Environ.* **2018**, *628*, 919-937. [CrossRef] [PubMed]

104. Sinnathamby, S.; Douglas-Mankin, K.R.; Craige, C. Field-scale calibration of crop-yield parameters in the Soil and Water Assessment Tool (SWAT). *Agric. Water Manag.* **2017**, *180*, 61-69. [CrossRef]

105. Gudmundsson, L.; Tallaksen, L.M.; Stahl, K.; Clark, D.B.; Dumont, E.; Hagemann, S.; Bertrand, N.; Gerten, D.; Heinke, J.; Hanasaki, N.; et al. Comparing Large-Scale Hydrological Model Simulations to Observed Runoff Percentiles in Europe. *J. Hydrometeor.* **2011**, *13*, 604–620. [CrossRef]

106. Kauffeldt, A.; Halldin, S.; Rodhe, A.; Xu, C.-Y.; Westerberg, I.K. Disinformative data in large-scale hydrological modelling. *Hydrol. Earth Syst. Sci.* **2013**, *17*, 2845–2857. [CrossRef]

107. Döll, P.; Kaspar, F.; Lehner, B. A global hydrological model for deriving water availability indicators: Model tuning and validation. *J. Hydrol.* **2003**, *270*, 105–134. [CrossRef]

<figure><img src="image_16.png" /></figure>

© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).