Independence of Future Changes of River Runoff in Europe from the Pathway to Global Warming

Abstract

The outcomes of the 2015 Paris Agreement triggered a number of climate impact assessments, such as for floods and droughts, to focus on future time frames corresponding to the years of reaching specific levels of global warming. Yet, the links between the timing of the warming levels and the corresponding greenhouse gas concentration pathways to reach them remain poorly understood. To address this gap, we compared projected changes of annual mean, extreme high, and extreme low river discharges in Europe at 1.5 °C and 2 °C under Representative Concentration Pathways RCP8.5 and RCP4.5 from an ensemble of regional climate model (RCM) simulations. The statistical significance of the difference between the two scenarios for both warming levels was then evaluated. The results show that in the majority of Europe (>95% of the surface area for the annual mean discharge, >98% for high and low extremes), the changes projected in the two pathways were statistically indistinguishable. These results suggest that in studies of changes at global warming levels, the projections of the two pathways can be merged into a single ensemble without major loss of information. With regard to the uncertainty of the unified ensemble, the findings show that the projected changes of annual mean, extreme high, and extreme low river discharge were statistically significant in large portions of Europe.

1. Introduction

River runoff is associated with different types of natural hazards. Extremely high flow regimes are associated with floods, a major natural hazard with considerable human and socio-economic implications [1]. Extreme low discharges are related to dry spells and droughts [2]. Changes in mean river discharge affect the long-term availability of water resources necessary for the sustainability of ecosystems and agricultural activities [3,4].

The past decades saw breakthroughs in atmosphere and precipitation modelling [5] and in our understanding of the complex dynamics of catchment areas [6] as well as in our observational [7] and computational capabilities, allowing an increasingly accurate investigation of large-scale, long-term issues and tendencies [1,8,9,10].

In view of global warming, the scientific community has put significant effort into studying how both mean and extreme river runoff values will change on the basis of climate projections, both at the global scale, such as in the Coupled Model Intercomparison Project 5 (CMIP5, [11]), and at the regional scale, such as in the Coordinated Regional Climate Downscaling Experiment (CORDEX, [12]), as well as to quantify the uncertainty in the projected changes [13]. Several past studies on future river runoff focused on specific periods in time, typically the short (2020s), medium (2050s), and long (2080s) term, in different representative concentration pathway (RCP) scenarios [14]. In such studies, the projected changes clearly depend on the considered pathway in, for example, studies related to river floods [15,16], droughts [17,18,19], and water resources [20,21]. It is worth mentioning the close relationship between the trends in river runoff and those of forcing variables such as precipitation and temperature. The projected changes in such variables have been extensively investigated [22,23] as well as their links with extreme and average river discharge [24,25].

In order to limit the impacts of climate change, the Paris Agreement [26] set the objective of limiting global warming to well below 2 °C and pursuing efforts to limit it to 1.5 °C compared to pre-industrial levels. This prompted scientists to study the changes in climate and hazards at global warming levels (GWL) of 1.5 and 2 °C as well as higher levels of warming to assess the potential consequences of not achieving the targets [27]. For example, this was assessed for floods [28,29,30,31], droughts [32,33], and for water resources [34].

The results of these studies rely on the hypothesis that the pathway to reach a certain greenhouse concentration and corresponding warming level plays a minor role in the change of the physical variables that define the hydrological hazard. The validity of such a hypothesis depends on the specific nature of the variables and on their constraint within a changing climate and should be verified in each case [35]. For example, for variables such as sea level rise, the pathway has been shown to play an important role [35,36] due to the large inertia of the sea masses. For average and extreme temperature and precipitation, there is no agreement on the degree of pathway dependence [35,37]. In particular, the physical mechanisms for precipitation are more complex than for temperature and involve many parameters, such as tropospheric moisture, aerosol concentration, and wind, and current knowledge of the underlying physics can only provide an order-of-magnitude estimation of the consequences of climate change [38]. Within these limitations, past contributions investigated the degree of dependency of temperature and precipitation from the pathway. Some of them suggest little effect of the pathway. For example, two studies analyzing the changes in temperature in Europe at different GWLs found between-scenario differences small compared with the climate model variability [39,40]. A similar study on precipitation found slightly more intense changes in mean precipitation in RCP4.5 than in RCP8.5 and attributed this to the difference in aerosol concentration between the two scenarios [41]. Other studies report significant differences among pathways, especially at the regional scale [42,43].

Whether or not the hypothesis of pathway-dependence to a global warming level is valid for hydrological variables that control riverine flood hazard and are connected with water availability and drought is yet unknown. This topic also relates with the investigation of the uncertainty associated with the emission pathway and of its comparison with other sources of variability. In this respect, this work contributes to the research on the quantification of epistemic uncertainty of modelling hydraulic variables and climate changes impacts and to the debate on the communication of the results [44,45].

In this study, we tested the pathway-dependence of future river runoff across a wide range of hydrological conditions. To this end, we produced a pan-European ensemble of river flow simulations with the Distributed Water Balance and Flood Simulation Model (LISFLOOD) [46,47] forced by a set of EURO-CORDEX [12] climate projections for the scenarios RCP8.5 and RCP4.5. The projected changes of extreme high, extreme low, and annual mean discharge between the baseline present climate and warming levels of 1.5 °C and 2.0 °C were thereby quantified. For each of these variables, the relevance of the between-pathway differences of projected change were evaluated with respect to the ensemble uncertainty. Finally, the two pathways were merged into a single ensemble, as has been previously suggested [48], and the statistical significance of the ensemble projected changes of extreme high, extreme low, and annual mean runoff was thereby quantified.

2. Materials and Methods

Projections of river discharge (Q) at a daily time step were obtained by running the hydrological model LISFLOOD [34,35] over the time horizon 1981–2100. The LISFLOOD represents all the relevant processes related to water dynamics in catchment areas including surface water balance (considering the effects due to the soil type, snow, frost, and vegetation), subsurface and ground water dynamics, river runoff and routing, and water demand for anthropogenic activities. For the sake of results reproducibility [49,50], it is worth mentioning that this model is now open source [51], and that version 2.8 was used in this study. The model domain covers the whole of Europe with a resolution of 5 km. The model parameters were calibrated with the European Flood Awareness System Meteo (EFAS-meteo) dataset from 1995 to 2015, optimizing the model output versus the streamflow data from 750 hydrological stations across Europe [28,52]. The forcing data of temperature, precipitation, radiative forcing, wind, and vapor pressure were provided by 11 bias-corrected [53,54] EURO-CORDEX regional climate models (RCMs) under scenarios RCP8.5 and RCP4.5 [12,53] (

Table 1. Ensemble members and related years when a global warming level (1.5 °C or 2 °C) was surpassed under Representative Concentration Pathways RCP8.5 and RCP4.5.

	RCP8.5, Year of Warming Level		RCP4.5, Year of Warming Level
CORDEX Model Name	1.5 °C	2 °C	1.5 °C	2 °C
CLMcom-CCLM4-8-17_CNRM-CERFACS-CNRM-CM5	2029	2044	2035	2057
CLMcom-CCLM4-8-17_ICHEC-EC-EARTH	2026	2041	2033	2056
CLMcom-CCLM4-8-17_MPI-M-MPI-ESM-LR	2028	2044	2034	2064
DMI-HIRHAM5-ICHEC-EC-EARTH	2028	2043	2032	2054
IPSL-INERIS-WRF331F	2021	2035	2023	2042
KNMI-RACMO22E-ICHEC-EC-EARTH	2026	2042	2032	2056
SMHI-RCA4_CNRM-CERFACS-CNRM-CM5	2029	2044	2035	2057
SMHI-RCA4_ICHEC-EC-EARTH	2026	2041	2033	2056
SMHI-RCA4_IPSL-IPSL-CM5A-MR	2021	2035	2023	2042
SMHI-RCA4_MOHC-HadGEM2-ES	2018	2030	2021	2037
SMHI-RCA4_MPI-M-MPI-ESM-LR	2028	2044	2034	2064

). The input data of potential evapotranspiration were estimated using the Evaporation Pre-Processor for LISFLOOD (LISVAP) model [55] with the Penman–Monteith parameterization [56]. In this study, we used version 0.3.2 of LISVAP which is now open source as well [57]. The simulations were carried out using the present estimations of anthropogenic land use, population, and anthropogenic water demand provided by the Land Use-based Integrated Sustainability Assessment of the Joint Research Centre (JRC LUISA) territorial modelling platform [58].

In order to understand how the emission pathway affects the projected changes in annual mean, and high and low extremes of Q, we analyzed for each model, scenario, and river pixel how these variables changed from the baseline (1981–2010) to the year of warming level (ywl), i.e., the year when a GWL was reached. For each model/scenario, ywl was estimated as the first year when the 30 year moving average of the GCM model’s global warming time series exceeded the GWL (Table 1).

The annual means discharge (Q_M) was estimated from the daily maps of Q. The extremes (high and low) of Q were computed by means of a non-stationary approach for extreme value analysis (EVA), the transformed-stationary EVA (tsEVA) [59]. This methodology uses discharge data over the whole time-horizon to fit the extreme value distribution (120 years in the present study), rather than data over the 30 year windows typically applied in stationary approaches. This implies that statistical fitting and extrapolation uncertainty is inherently lower which is especially relevant for extreme events with return periods that go beyond 30 years such as the 100 year flood. Furthermore, non-stationary techniques are better able to capture the changing statistics at different warming levels, especially in RCP8.5, where the warming continues after exceeding 2.0 °C. The authors of tsEVA released an open source MATLAB toolbox [59] which was employed in this study.

For high extremes, tsEVA was set-up to fit the peaks of Q beyond a moving 98.5 percentile in 1981–2100 with a non-stationary generalized pareto distribution (GPD) thus estimating the time-varying 100 year return level Q_H100 and its uncertainty. For each pixel, the input of tsEVA was the daily series of Q, and the algorithm was set-up in order to consider peaks at a distance of at least 30 days from one another.

The analysis of the low extremes was, to some extent, complicated by the fact that Q is low-bounded to 0 and that a 0 runoff condition is met regularly in many points. This can lead to an unreliable fit of the extremes (e.g., GEV distributions with shape parameter < −0.5) and to an estimation of the return values too close/equal to zero which would result in inaccurate projections of relative changes. To solve this problem:

(a)

The daily time series of Q were smoothed with a monthly running mean before fitting an extreme value distribution, guaranteeing that medium-term lows were considered in the fit;

(b)

The generalized extreme value (GEV) distribution was used instead of the GPD, due to the excessive proximity of the GPD threshold to 0 in many pixels;

(c)

The pixels with a fitted GEV shape parameter < −0.5 were excluded from the analysis;

(d)

The time-varying 15 year low return level Q_L15 was preferred to the 100 year one, due to the excessive proximity of the latter to 0 in many points.

For Q_H100, Q_L15, and annual mean discharge Q_M, the relative projected change ΔQ between the baseline period (from 1981 to 2010) and the thirty years centered in ywl was computed for each scenario, RCM, and warming level:

$$\mathsf{\Delta }{Q}_x=\frac{Q_{x,wl}-Q_{x,bsl}}{Q_{x,bsl}}$$

(1)

where the suffix “x” indicates one of the three considered quantities (Q_H100, Q_M, and Q_L15), and bsl and wl indicate the values at the baseline and at the warming level, respectively.

For each GWL and analyzed variable, the two scenarios were joined into a single 22 model ensemble, following an approach suggested by previous studies [48]. This allows for a better quantification of climate uncertainty that also accounts for the uncertainty associated with the emission pathway. The variance of the ensemble was studied, pixel by pixel, by means of a one-way, two-group analysis of variance (ANOVA, e.g., [60,61]). This approach provides information on the significance of the difference between the two groups, similar to a Welch’s t-test [62], but ANOVA was preferred, as it comes with a decomposition of the variance of the projected change, σ², into its components due, respectively, to the inner (or intramodal) variability of the two pathways and to the difference between the pathways:

$${\sigma }^2={\sigma }_{within}^2+ {\sigma }_{between}^2$$

(2)

The importance of the between-scenario difference with respect to σ_within was evaluated by means of an F-test: the pathways are considered as significantly different where the p-value ≤ 0.05. The ANOVA analysis was carried out using the python library scipy.stats version 1.4.1 [63].

Finally, the statistical significance of the ensemble projected change of Q_H100, Q_L15, and Q_M was studied, classifying the change as significant when |ΔQ|>σ as proposed in a previous study on floods [16].

3. Results and Discussion

The ensemble changes projected for Q_H100, Q_M, and Q_L15 at 1.5 °C and 2.0 °C were similar in the two pathways with projected changes generally more intense at 2 °C than at 1.5 °C. Changes of Q_H100 were in the range −15% to 32% (values at the 99th and the 1st percentiles) at 2 °C with respect to the baseline (Figure S1). The ΔQ_M values were between −19% and 36%. Q_L15 was the variable that showed the strongest changes, ranging between <−50% to >100% (Figure 1), in line with the fact that small quantities are subject to larger relative variations.

The results of the ANOVA show that the difference between the two pathways is generally small compared with the inter-model variability for all of the three variables. On average, the component of the variance σ²_between due to the differences between scenarios explains ≤4% of the total variance for all of the three considered variables (≤3% for high and low extremes) with slightly lower percentages at 1.5 °C than at 2.0 °C (

Table 2. Components of the total standard deviation for the considered ensembles expressed as percentages of explained variance.

		% Explained Variance		% Points Where p-Value ≤ 0.05
Variable	Warming Level	σ2within	σ2between
QH100	1.5 °C	98.6%	1.4%	0.1%
	2.0 °C	97.0%	3.0%	2.0%
QM	1.5 °C	98.2%	1.8%	0.2%
	2.0 °C	96.0%	4.0%	4.4%
QL15	1.5 °C	97.2%	2.8%	0.7%
	2.0 °C	97.1%	2.9%	1.2%

). For all of the three variables and for both warming levels, the intra-scenario standard deviation σ_within (Figure 2c,d, Figure 3c,d, and Figure 4c,d) generally had the same order of magnitude of the projected change (Figure 2a,b, Figure 3a,b, Figure 4a,b), while σ_between was much smaller (Figure 2e,f, Figure 3e,f, Figure 4e,f). The ANOVA F-test showed that the differences among scenarios were generally statistically non-significant: at 1.5 °C, the p-value was larger than 0.05 in >99% of the points for the three considered variables (99.9% for Q_H100, 99.8% for Q_M, and 99.3% for Q_L15). At 2.0 °C, these percentages decreased, but remained >95% (98% for Q_H100, 95.6% for Q_M, and 98.2% for Q_L15; Table 2). We can further argue that in this ensemble, the inter-model variability was solely related to the CORDEX forcing data, while the setup of LISFLOOD was constant across all the simulations. In other words, we are not taking into account the epistemic uncertainty arising from our lack of knowledge of catchment dynamics [64,65]. Adopting different setups or even using different hydraulic models may result in increased inter-model variability.

The differences between the two scenarios, though generally small, point toward more intense changes in RCP8.5 for all of the three variables, especially at 2 °C. This result is driven by similar patterns of change in precipitations (not shown). However, global studies on precipitation report a more intense increase in precipitation in scenario RCP4.5, explaining it by the fact that the concentration of aerosol is projected to decrease approximately with the same annual rate in the two scenarios. Consequently, when the warming level is reached, the sky is clearer of aerosol in scenario RCP4.5 [41], and this comes with a weaker shortwave radiation absorption and reflection by the atmosphere, an increased amount of shortwave radiation reaching the land surface, an increased top atmosphere radiative cooling and subsequent more intense precipitations [66]. The relevance of aerosol concentration for explaining between-pathway differences in precipitations is remarked upon also by other authors [43,67,68]. The results of these studies are global and, therefore, do not contrast with our finding of slightly more intense precipitations in scenario RCP8.5 over Europe. The projected decreased rate in aerosol concentration is not uniform across the world and across scenarios, and the differences we found in Europe might be associated with a faster depletion of black carbon and sulfur aerosols in RCP8.5 by the end of the century [69].

Overall, these results support the assumption for merging the two pathways into a single 22 model ensemble, rather than studying separately the two scenarios, for a more comprehensive estimation of extent and uncertainty of the projected changes at the two warming levels.

In the majority of Europe, the change projected by the merged ensemble was positive for the three variables (i.e., higher extremes and annual means). At 1.5 °C ΔQ_H100, ΔQ_M, and ΔQ_L15 were positive in 83%, 94%, and 77% of the pixels, respectively. At 2 °C, these figures slightly decreased to 80%, 89%, and 73%. Accordingly, the spatial median of the projected change was positive for all of the three quantities: at 1.5 °C the median changes were ~8%, ~10%, and ~15% for Q_H100, Q_M, and Q_L15, respectively. At 2 °C these figures increase up to ~10%, ~14%, and ~38%.

The spatial variability of the changes differs for the three variables. For Q_H100, the changes were positive in central Europe and in the majority of southern Europe, with the exception of Andalusia, part of the southern Balkans, and Sicily, and which were negative in the majority of Scandinavia (Figure 2a,b). For Q_M, the changes were positive in central and northern Europe, turning negative in the Iberian Peninsula, and to a lesser extent in the southern part of the Italian Peninsula and Greece (Figure 3a,b). For Q_L15, the differences between northern and southern Europe were the strongest, with an intense increase in extreme lows in north-eastern Europe and an intense decrease in the south-west, the two regions separated by a slow gradient (Figure 4a,b).

The statistical significance of the projected changes was different for the three variables. ΔQ_H100 was significant in 19% and 24% of the pixels at 1.5 °C and 2.0 °C, respectively (Figure 2a,b). Significant positive changes of ΔQ_H100 were estimated for central Europe, covering the larger part of France, Germany, the British Isles, Benelux, Denmark, the northern/central Danube basin, western Poland, and the Po River valley, while the intense changes projected in eastern Europe and in the Dniepr basin were subject to strong uncertainties. Significant negative changes in ΔQ_H100 were found in Iceland and in part of Scandinavia.

The ΔQ_M was significant in 50% and 54% of the pixels at 1.5 °C and 2.0 °C, respectively (Figure 3a,b). The area of positive significant changes included the whole central–north–eastern Europe and Iceland, excluding Romania and of part of Ukraine. From 1.5 °C to 2.0 °C, a decrease in the significance of changes in the British Isles and in northern France can be observed, while the positive changes in eastern Europe and the negative changes in Andalusia become more significant.

The ΔQ_L15 was significant in 33% and 41% of the pixels at 1.5 °C and 2.0 °C, respectively (Figure 4a,b). A significant increase in Q_L15 can be observed in north-eastern Europe and in Iceland, while a significant decrease is shown in the Iberian Peninsula, parts of France, Italy, and the Balkans with changes more significant at 2.0 °C than at 1.5 °C. The projected changes of annual means Q_M were significant for a larger area than the extremes Q_H100 and Q_L15. A possible explanation for this is that the trend of Q_M follows that of the annual mean precipitation projected by the pathways. While Q_H100 and Q_L15 were associated with extreme meteorological conditions that are inherently more challenging to identify [70].

Overall, the results on the projected changes and on their statistical significance confirm previous results. A study from 2012, carried out with an ensemble of 12 models in the emission scenario SRES A1B, reports similar patterns of change for annual mean runoff, although with stronger negative changes in southern Europe [15]. The main difference in the projection of the 100 year return levels was a reversal of sign in the projected change in eastern Europe. Other authors, using a smaller CORDEX ensemble in scenario RCP8.5, found patterns of change similar to ours for both annual means and high extremes [16]. A further study, carried out with a different CORDEX ensemble with three RCP scenarios, identified patterns of change of low flow similar to the ones found in this study for QL15 [32]. Previously, a study on projected river flow regimes in 2050, carried out with the model WaterGAP3 forced by three GCMs found patterns similar to the ones described in this contribution [71]. Other studies and reviews reporting qualitatively similar results can be mentioned [25,29,72,73,74,75].

Merging the two pathways into a single ensemble comes with a two-fold advantage with respect to the separate treatment of the two scenarios. On the one hand, it improves the estimation of the statistical significance of the projected change by increasing its size from 11 to 22 and by better taking into account the pathway-related uncertainty (the emission pathways are set ex ante as a hypothesis for the CMIP experiment and are generally not considered as a source of uncertainty). On the other hand, a multi-pathway ensemble can simplify the discussion of the projected changes by removing from the analysis the dependency from the emission pathway and making the results clearer and more understandable by a non-scientific public.

On the other hand, a limitation of a multi-scenario ensemble is that each emission pathway reaches similar projected changes at different time frames. But vulnerability and exposure to climate hazard changes over time [76]. Therefore, the application of multi-scenario ensembles to studies that include evolving socio-economic scenarios (e.g., to model the effects of climate extremes) or focusing on mitigation or adaptation would require further assumptions for the quantification of the time available for acting.

4. Conclusions

In this study, we examined the high and low extremes and the annual means of river discharge in a changing climate at warming levels of 1.5 °C and 2.0 °C under greenhouse emission scenarios RCP8.5 and RCP4.5. These three variables were chosen to provide full coverage of different regimes that have impacts on human activities, ecosystems, hazard, and risk. The variance of the projected changes was studied with the ANOVA methodology to investigate the statistical significance of the differences between the two scenarios. Our results indicated that the projected changes were approximately a function of the warming level, as the between-pathway differences were generally much smaller than the within-pathway variability. The projected changes and their uncertainties were thereby studied on the joint two-pathway, 22 model ensemble for both warming levels. The changes in high extreme runoff were found to be generally positive in central Europe and negative in Scandinavia and in Southern Europe. The other two variables exhibited changes mostly positive in Northern Europe and negative in Southern Europe. For all three variables, the changes were generally more intense and statistically significant at 2.0 °C than at 1.5 °C. These results roughly agree with previous single-scenario studies.

A multi-scenario ensemble was thus found to be a viable way to provide additional information on future changes at GWLs compared to analyzing separately the different scenarios. It allows a more reliable evaluation of the total climate uncertainty: the future emission pathway is itself unknown, not only due to the complexity of climate dynamics and the uncertainty therein, but also due to the impossibility of predicting the commitment of the international community in mitigation efforts or the effects of such efforts. A multi-scenario ensemble approach to portray impacts at GWLs is a way to quantify the variability associated with the emission pathway alongside the within-pathway climate variability.

Last but not least, a multi-scenario ensemble can simplify the way future climate change projections and effects of global warming levels are communicated to the public. This is relevant for climate policy [45]. Shortening the discussion on the effects of different emission pathways can contribute to providing a clearer, more compact, and summarized account of the results. It thus facilitates to vehicle the information to stakeholders and a non-scientific public and, therefore, favors the dissemination of knowledge.