Compensatory climate effects link trends in global runoff to rising atmospheric CO₂ concentration

Observed monthly river discharge records, for years 1960–1999 inclusive, are used for JULES-CN model evaluation. Discharge measurements are from the farthest downstream stations, and for 42 large basins around the world (figure S1 is available online at stacks.iop.org/ERL/14/124075/mmedia). These are obtained from three sources: (i) Global River Discharge Centre (http://bafg.de/GRDC); (ii) China Statistical Yearbook (Bureau 2000); and (iii) river discharge archive by Dai et al (2009). The two main selection criteria for catchments are they need to have >80% available data for the study period, and catchment areas are >100 000 km² to match the relatively large spatial resolution (gridbox of 1.25° lat × 1.875° lon, ∼60 000 km²) used in the JULES-CN simulations of this study.

The Joint UK Land Environment Simulator (JULES; Best et al 2011, Clark et al 2011) is the land surface scheme of the UK Met Office Earth System Model. JULES can also be operated 'offline', as here, and forced with known surface meteorological conditions for the contemporary period. JULES evolved from the Met Office Surface Exchange Scheme (MOSES; Cox et al 1998, 1999), combined with a dynamic vegetation module called the Top-down Representation of Interactive Foliage and Flora Including Dynamics (TRIFFID; Essery et al 2003, Clark et al 2011). The current JULES version contains a sophisticated representation of canopy radiation interception (Mercado et al 2007) that defines explicitly the diffuse and direct components of the photosynthetically active radiation (Mercado et al 2009). It also includes a representation of the response of vegetation to tropospheric ozone deposition (Sitch et al 2007).

Notably, JULES-CN has implemented the nitrogen cycle into the dynamic global vegetation and hydrology components of the JULES model (Wiltshire et al 2019). In particular, processes of soil nitrogen dynamics, litter production, plant uptake, nitrogen allocation, response of photosynthesis and maintenance respiration to varying nitrogen concentrations in plant organs and inorganic nitrogen are now included in the model. This addition improves the simulation of vegetation distribution, as it responds to climate, CO₂ and nutrient availability. This new JULES framework also allows an assessment of the impacts of vegetation carbon–nitrogen interactions on hydrological processes, including those driven by changing atmospheric nitrogen deposition. In brief, the nitrogen module acts to reduce Gross Primary Productivity (GPP) in regions where insufficient inorganic nitrogen exists to meet demand. In regions of nitrogen limitation, it is expected that the vegetation biomass and leaf area will be less stimulated under raised atmospheric CO₂ concentrations, compared to the model version without nitrogen representation (i.e. JULES-C). The effect of nitrogen is to generally cause altered vegetation structure, as opposed to changes to stomatal conductance. Inorganic sources of nitrogen include deposition, fixation and soil net mineralization, and in some circumstances, these offset nitrogen limitation.

We run two sets of simulations: one using the JULES model with carbon–nitrogen interactions (JULES-CN) and one without carbon–nitrogen interactions (JULES-C). Such simulations, named 'factorial', isolate the individual contributions of elevated CO₂ concentration, climate change and LUC predicted by JULES-C and JULES-CN models (table 1). Here, JULES-C and JULES-CN simulate historical land surface conditions, driven by a 99 year (1901–1999) observation-based meteorological forcing dataset, CRU-NCEP v4 (Le Quéré et al 2012). In addition, annual atmospheric CO₂ concentration data is from Keeling and Whorf (2005), LUC data from the HYDE dataset (Klein Goldewijk 2011) and nitrogen deposition from ACCMIP (Lamarque et al 2013). All simulations are 'spun up' to equilibrium (i.e. at 1901) under environmental conditions by a repeated 20 year climate forcing data (1901–1920) and fixed pre-industrial values of CO₂, land-use and nitrogen deposition.

Comprehensive global observations of diffuse and direct partitioning of downward shortwave radiation and ozone concentrations are unavailable. The simulations in table 1 use the fixed diffuse fraction of 0.4 and the pre-industrial O₃ concentration. To assess O₃ and aerosol radiative effects on runoff changes, we add an extra second set of two JULES-CN simulations (table 2). We use HadGEM2-based estimates of diffuse fraction and O₃ concentration, and combine with the same estimates of surface meteorological conditions that are used for the simulations in table 1 as the forcing. The radiative effects of atmospheric aerosols adjust the balance of downward shortwave between direct and diffuse levels. Here, the aerosol forcing data, similar to Mercado et al (2009), includes tropospheric aerosol species: black carbon, sulphate, mineral dust, sea salt and biomass burning. Distributions of aerosol optical depths for each species were simulated by the atmospheric model in HadGEM2-A, following Bellouin et al (2007). Radiative transfer calculations are used to estimate diffuse fraction at the land surface, based on these aerosol distributions, which surface level is input to JULES (Mercado et al 2009). The tropospheric ozone forcing data, similar to Sitch et al (2007), is generated by 'time-slice' simulations with a tropospheric chemistry model STOCHEM (Sanderson et al 2003). Values for the intermediate years are generated by linear interpolation.

To illustrate how the dominant drivers vary with spatial scale, we follow the method used by Jung et al (2017) and define the relative magnitudes of environmental components ('COMP') that force runoff trends. This metric is defined as the trend of runoff forced by an individual component, divided by the trends of runoff forced by all factors (trend_ALL): ${D}^{{\rm{C}}{\rm{O}}{\rm{M}}{\rm{P}}}=\left|{\rm{t}}{\rm{r}}{\rm{e}}{\rm{n}}{{\rm{d}}}_{{\rm{C}}{\rm{O}}{\rm{M}}{\rm{P}}}\right|/\left|{\rm{t}}{\rm{r}}{\rm{e}}{\rm{n}}{{\rm{d}}}_{{\rm{A}}{\rm{L}}{\rm{L}}}\right|.$ This index is calculated for spatial windows of 1 × 1, 2 × 2, 4 × 4, 8 × 8, 16 × 16, 28 × 32, 56 × 48, and all the grid cells over the globe. The magnitude of spatial scales (in km) defined by the longitudinal distance at 45^oN is also shown in figure 7.

We first evaluate the JULES model performance in reproducing historical river discharge for period 1960–1999 (figure 1). The JULES-CN model performs well in capturing the long-term averages of annual discharge (figure 1(a)) and the inter-annual variability of annual discharge (figure 1(b)) for the 42 major basins, and for the period 1960–1999, when the effects of all six forcings are included.

Zoom In Zoom Out Reset image size

The JULES-CN model performs well for many basins e.g. Congo, Mississippi and Yangtze. However, it has a limited ability to reproduce observed runoff increases in a few basins, e.g. some Eurasian Arctic rivers, the Amazon and Pakistan basins. Increased precipitation plays a major role in observed Eurasian river discharge increases, but poor performance there could be partly related to the uncertainty in regional precipitation. For instance, the NCEP precipitation data shows smaller increases than other precipitation datasets in high latitudes (a point also noted in Pavelsky and Laurence 2006). In addition, these Eurasian Arctic rivers are highly affected by permafrost (figure S2), which as yet is not fully included in JULES-CN. The bias in trends for Asian basins may reflect poor representation of direct human intervention (notably dams, irrigation) on runoff (Yang et al 2018). Raised levels of irrigation and regulation enhance surface evaporation and possibly reduce runoff in intensively cultivated areas (hence our findings support those of Tang et al 2007, Haddeland et al 2014). In the Amazon basin, the model-data discrepancy might be partly related to an inadequate representation of the response to land use, notably deforestation.

There are strong similarities between figures 2(b)–(d), providing initial evidence that AER and O3 effects are likely relatively small drivers.

Zoom In Zoom Out Reset image size

Annual global runoff averages are calculated from the 42 large river basins (figure S1), for comparison with observations. The time-evolving, observation-based global average of runoff has a small (but statistically significant) positive trend of +0.22 ± 0.20 mm yr⁻² (+0.07 ± 0.06% yr⁻¹) during the period of 1960–1999 (figure 3(a)). The JULES-CN simulation CN_{CO2+CLIM+LUC+NDE} provides a comparable estimate of the global runoff trend of +0.30 ± 0.27 mm yr⁻².

Zoom In Zoom Out Reset image size

Further analyses use the factorial simulations to quantify the contributions of six environmental factors to the global averages of modelled runoff linear trends for the period of 1960–1999 (figure 3). Values are also calculated from the 42 large study basins, and are first presented as the accumulation of each factor (yellow bars, figure 3), before showing the isolated effects of the individual drivers. Contributions from changing climate are small for the large-scale averages (just +0.01 ± 0.30 mm yr⁻²). Rising atmospheric CO₂ concentration is the most important contributor to the global river runoff trend, with a significant increase of +0.18 ± 0.006 mm yr⁻². Important contributions are also LUC (increase in total runoff of +0.05 ± 0.007 mm yr⁻²) and carbon–nitrogen interactions and nitrogen deposition (decrease of –0.03 ± 0.007 mm yr⁻²). The radiative effects of atmospheric aerosols adjusting the balance of downward shortwave between direct and diffuse levels, and ozone concentration changes, both exert a slightly negative influence on global runoff changes. Thus, there are no significant differences among the globally-averages of runoff change between CN_{CO2+CLIM+LUC+CN&NDE}, CN_{CO2+CLIM+LUC+CN&NDE+AER} and CN_{CO2+CLIM+LUC+CN&NDE+O3} simulations, supporting the noted small effects on global runoff trends of aerosol and O3 forcings (figures 2(b) versus (c) and (d)).

Figure 3 shows at global scale, almost all drivers have relatively little effect on runoff trends, except for increasing atmospheric CO₂. However, the findings may be different regionally, as suggested by figure 2(a), which shows not only spatial heterogeneity but also large differences in sign. We return to consider geographical differences, and find that at the local level, JULES-CN simulates a significant increase in runoff during the late twentieth century over eastern United States, northern Europe, some areas in China and a large area in South America. Decreased runoff is simulated over some areas of Canada, central Africa and west Amazon (figure 4(b)). These changes compare well to data (figures 2(a) and 4(a)).

Zoom In Zoom Out Reset image size

For local runoff trends, multiple factors show much more diversity in their relative contributions to runoff changes. The driving factor of the largest magnitude and for each grid cell, as calculated by JULES-CN factorial simulations, is presented in figure 4(c). Figure 4(c) shows that unlike the mean dominant contributors to the global runoff trends, climate change is the factor responsible for the absolute largest trends of runoff. Climate as the dominant driver is for over 82% of global land area (excluding Antarctica), whereas rising atmospheric CO₂ concentration dominates runoff increases for <5% of global land area. For locations where climate change dominates, climate change impact on runoff is positive over 38% of global land area, and negative for the remaining 44%. For the remaining runoff increasing areas where climate is not the dominant driver, the main driving factor is either LUC or tropospheric ozone. This is especially notable for some regions in China and India (figure 4(c)). We also identify the factor with the second-largest absolute value of runoff trend, and observe this has strong geographical heterogeneity (figure 4(d)). Elevated CO₂ concentration, via its physiological effects, contributes to an upward trend in the runoff, and is the second-largest effect over 35% of global land area. In contrast, carbon–nitrogen interactions and nitrogen deposition induced the second-largest decreasing trends in runoff over the nitrogen-limited regions, e.g. Australia and the boreal regions at the high latitudes. LUC contributes to the second-largest runoff increase in west Amazon.

Actual trend values for different drivers are presented in figure 5. Trends in simulated runoff due to climate change, over the period 1960–1999, are shown in figure 5(a). The climate change signal (figure 5(a)) is predominantly a function of increasing atmospheric CO₂ concentration, which is the dominant driver when expressed via contribution to radiative forcing (IPCC 2013). Other smaller drivers are non-CO₂ greenhouse gases, volcano effects, solar fluctuations and atmospheric aerosols. Hence, figure 5(a) is regarded as broadly, an indirect CO₂ effect. Besides temperature increases, rainfall intensity increases in some areas, wind speed changes affect near-surface turbulent exchange and surface conditions, which in turn can adjust evapotranspiration. This also causes an indirect effect, where the response to changed surface conditions alters land-atmosphere CO₂ exchange, and ultimately vegetation dynamics (McVicar et al 2012, Dourte et al 2015) (figure S3). However, the climate-induced runoff trends are correlated significantly with the trends of precipitation over 63% of global land area (figure S3), suggesting precipitation is the dominant local climate driver.

Zoom In Zoom Out Reset image size

The net physiological and structural effects of vegetation response and impact on runoff, caused by rising atmospheric CO₂ concentration, is shown in figure 5(b). Higher CO₂ concentrations generally enhance runoff in humid regions, whereas it slightly reduces runoff in arid regions. For the wet regions, JULES-CN simulations suggest that runoff increases are related to mostly reductions in transpiration caused by CO₂-induced stomatal closure (figure S4(b) and text S1) reducing transpiration (figure S4(d)). In water-limited environments, the elevated CO₂ concentration reduces stomatal conductance instantaneously, through stomatal closure, to alleviate plant dryness stress. Howver, over longer time periods (e.g. decades), elevated CO₂ concentration leads to increased vegetation WUE, and which triggers increases in LAI and/or extended growing seasons. Such increases in LAI can enhance transpiration, through increased stomata density and number (figure S4(c)), which offset the retention effects of stomata closure (figure S4(d)). Additionally, higher LAI can raise canopy interception and related evaporation increases (figure S4(f)). The balance of all of these effects in JULES is to lead to a slight decrease in runoff.

Figure 5(c) shows the impact of LUC, leading to runoff increases over the wet tropical regions (especially in the Amazonian region) and decreases in Southeast China, Central Asia, the eastern USA. The net losses of forest area (due to deforestation and/or fire; Klein Goldewijk et al 2011) and cropland expansion during the past four decades play multiple important roles. In the Amazonian region, simulations with JULES-CN model suggest that when forest is replaced by cropland (figures 6(a)–(c)), ET decreases, resulting in overall runoff increases, since compared to cropland (or grassland), forest has higher transpiration rates (Mahmood and Hubbard 2003). In contrast, the replacement of croplands in Southeast China is more productive than the replaced forest (figure 6(d)), leading to an increase in ET and a decrease in runoff. Also, runoff in Central Asia has decreased due to the local expansion of agriculture. Interestingly, in the eastern USA, a vegetation shift from C4 grasslands towards C3 grasslands has likely occurred (figures S5(c) and (d)). C3 grasslands have low WUE, which implies their ability to maintain photosynthesis depends on using more water, thus leading to decreased runoff.

Zoom In Zoom Out Reset image size

Nitrogen availability limits plant growth, particularly over the temperate and boreal regions (LeBauer and Treseder 2008), potentially offset by raised nitrogen deposition. Explicitly modelling carbon–nitrogen interactions, JULES-CN simulations suggest that increased nitrogen deposition in Siberia, Middle East and southwest China results in increasing GPP, LAI and evapotranspiration, thereby decreasing runoff (figure 5(d)). However, in the Eastern Europe and North America, nitrogen limitation effects dominate, leading to decreasing leaf photosynthesis and transpiration, and thereby increasing runoff (figure S6).

Figure 5(e) shows that the impact of increases in aerosols in densely-populated and in biomass burning regions (figures S7(a) and (b)) during 1960–1999. This includes South America, central Africa, Australia, India, southern and western Asia. In contrast, most of south-west USA and western Europe experience aerosol decreases due to Clean Air policies. The upward trend of runoff in western Europe (except the Iberian Peninsula) and south-west US demonstrate the effect of diffuse radiation fraction (Mercado et al 2009 for details of diffuse fraction changes) decreases are larger than the opposing influence of total radiation increase. Raised atmospheric aerosol levels lead to diffuse fraction increases in South America, central Africa, Australia, India, southern and western Asia, which again dominant direct radiation changes. For these locations, this leads to increased transpiration and decreased runoff.

Finally, simulated changes in tropospheric ozone concentration over the study period leads to increased simulated runoff (figure 5(f)) with changes being statistically significantly over most of South America, Africa and southern China. Despite statistical significance, the magnitude of runoff trends due to ozone changes is relatively small.

Investigated further is our finding that at the local scale, climate change is the dominant driver of runoff trends, whereas when integrated globally, the dominant factors are elevated CO₂ concentration and LUC. We explore whether compensatory effects of climate change when scaling from local to global scales explain this paradox. The calculated relative magnitude of climate-induced runoff trend (D^CLIM) decreases with increasing spatial aggregation, while the relative magnitude of CO₂-induced runoff trend (D^CO2) increases (figure 7). The decreases in D^CLIM is due to a compensation of positive and negative trends of climate-induced runoff between different grid cells (as shown in figure 5(a)), whereas the CO₂-induced runoff trends are generally universal and positive (figure 5(b)). This strong dependence of main driver on scale has similarities to findings by Jung et al (2017), but with their analysis in terms of land-atmosphere CO₂ exchanges.

Zoom In Zoom Out Reset image size

As outlined above, rising atmospheric CO₂ concentration influences plants and associated water cycle in two contrasting ways. First, the CO₂ fertilization effect changes vegetation structure, stimulates photosynthesis and raises the biomass of C3 plants. This first effect has three consequences: (i) it increases transpiration by raising the number of stomata due to increased LAI; (ii) it increases canopy interception loss due to raised LAI and (iii) soil evaporation decreases due to reduced available energy at the surface. Such changes reduce local runoff. Second, there is the direct physiological effect, causing stomatal closure, (although higher WUE), leading to transpiration decrease and runoff increase. JULES-CN calculate the balance between these two-opposing effect is such that runoff in drylands has decreased, but in wet regions it has increased. Gerten et al (2008) performed a similar analysis with the LPJmL model, which is a C-only model. They found rising atmospheric CO₂ levels decreased runoff in some drylands and increased runoff in temperate and boreal regions, which is consistent with JULES-CN. However, LPJmL simulated non-significant trends in runoff over the tropics, which are smaller changes than we predict. Using the ORCHIDEE model (C-only model), Piao et al (2007) found the net effect of CO₂ on runoff to be overwhelmed by CO₂ fertilization effects (i.e. rather than direct stomata closure effects), and so negative across tropical wet and temperate regions, and thus different to our findings. In addition, using 13 land surface models, Mao et al (2015) found at elevated CO₂, for most areas and especially these regions covered by tropical broadleaf evergreen trees and high latitude shrubs, showed decreasing trends in evapotranspiration, whereas dry areas with sparse vegetation showed increasing evaporation. This result is consistent with our findings, since runoff often changes in opposite directions to evapotranspiration (since runoff is equal to precipitation minus evapotranspiration and precipitation shows relatively little change for altered CO₂).

In addition to the uncertainty associated with direct and indirect CO₂ influences, until recently land surface models have lacked representation of the terrestrial carbon–nitrogen cycle. Terrestrial carbon uptake is modelled as limited by the availability of nitrogen in JULES-CN simulation. Specifically, an ecosystem becomes limited when insufficient nitrogen is available for plants to allocate net photosynthetic to growth, in which case photosynthesis is 'downregulated' to match the available inorganic nitrogen. Nitrogen does not directly affect photosynthetic capacity through leaf nitrogen concentrations, but instead acts indirectly by controlling the biomass and LAI. An increase in nitrogen limitation would generally act to suppress the evaporative fraction and vice versa for a decrease in limitation. In other words, CO₂ fertilization is known to be constrained by nitrogen availability (Felzer et al 2009, Norby et al 2010, Zaehle 2013), and therefore models without a coupled C and N cycle may overestimate fertilization-related runoff decreases (Hungate et al 2003, Zaehle and Dalmonech 2011). We compare projections of runoff from JULES-CN (i.e. with carbon–nitrogen interaction processes) and JULES-C model (without carbon–nitrogen interaction processes). We find JULES-C projects a weaker increase in CO₂-induced runoff trends over most regions of the world, especially for boreal forest regions, compared to the JULES-CN model (figure 8(b)). This is expected, as JULES-C model projects larger CO₂ fertilization-induced runoff decreases, which offsets more runoff increases due to raised stomatal closure. Comparison of figures 8(a) and (b) shows the relative importance of the inclusion of the N-cycle on runoff. We note that our results are specific to the JULES model. We hope our analysis will act as an incentive for other DGVM groups, to determine the impact of new geochemical cycle modelling (and especially the nitrogen dynamics) on land impacts of concern, and in particular runoff.

Zoom In Zoom Out Reset image size

Our projections do have caveats. Our simulations are in 'offline' mode (forced by reanalysis data), hence lacking any land-atmosphere feedbacks in response to changes in land surface configuration (Smith et al 2016); this might modulate future predictions of runoff. Uncertainty remains in representation of several eco-hydrological processes, and including stomatal conductance-photosynthesis coupling and the transpiration response of plants (De Kauwe et al 2013, Swann et al 2016). Land models are known to exhibit strong differences and biases in their aggregation from leaf to full canopy fluxes (Lian et al 2018).