Global runo ff over 1993 – 2009 estimated from coupled land-ocean-atmosphere water budgets and its relation with climate variability

Introduction Conclusions References


Abstract
Whether the global runoff (or freshwater discharge from land to the ocean) is currently increasing and the global water cycle is intensifying is still a controversial issue.Here we compute land-atmosphere and ocean-atmosphere water budgets and derive two independent estimates of the global runoff over the period 1993-2009.Water storage variations in the land, ocean and atmosphere reservoirs are estimated from different types of datasets: atmospheric reanalyses, land surface models, satellite altimetry and in situ ocean temperature data (the difference between altimetry based global mean sea level and ocean thermal expansion providing an estimate of the ocean mass component).Results for the global runoff from the two methods show a very good correlation between both estimates.More importantly, no significant trend is observed over the whole period.Besides, the global runoff appears to be clearly impacted by large-scale climate phenomena such as major ENSO events.To infer this, we compute the zonal runoff over four latitudinal bands and set up for each band a new index (Combined Runoff Index) obtained by optimization of linear combinations of various climate indices.Results show that, in particular, the intertropical and northern mid-latitude runoffs are mainly driven by ENSO and the Atlantic Multidecadal Oscillation (AMO) with opposite behavior.Indeed, the zonal runoff in the intertropical zone decreases during major El Ni ño events whereas it increases in the northern mid-latitudes, suggesting that water masses over land are shifted northward/southward during El Ni ño/La Ni ña.
In addition to this study, we propose an innovative method to estimate the global ocean thermal expansion.The method is based on the assumption that the difference between both runoff estimates is mainly due the thermal expansion term not accounted for in the estimation of the ocean mass.Comparison of our reconstructed thermal expansion with two existing datasets shows the relevance of this new method.

Introduction
Continental waters are continuously exchanged with atmosphere and oceans through vertical and horizontal mass fluxes (precipitation, evaporation, transpiration of the vegetation, surface runoff and underground flow).Freshwater discharge from land to ocean (or global runoff) is a key component of the global water cycle.Its variability reflects the continental hydrological dynamics and is then impacted by climate change (e.g.intensification of precipitation over land) and anthropogenic activities (reservoirs, land use changes, irrigation, groundwater pumping).As noted by many authors (e.g.Labat et al., 2004;Huntington, 2006;Gerten et al., 2008;Dai et al., 2009), global runoff may be seen as an indicator of the intensification of the hydrological cycle.Besides, since global runoff represents an integrated response to continental hydrological dynamics, it has also been used to detect the impact of anthropogenic activities (e.g.Gleick, 2003;Nilsson et al., 2005;Milliman et al., 2008).
Numerous studies have focused on the characteristics of global runoff, in terms of long term mean, trends, spatial distribution and interannual variability.The most basic way to estimate global runoff is the use of gauged based measurements of river discharge at the outlet of the world major basins.Data bases such as the Global Runoff Data Centre (GRDC) provide such data for a large number of gauge stations worldwide.Although this approach gives the most direct estimation of global runoff, it remains limited by some important drawbacks (see e.g.Legates et al., 2005;Peel and McMahon, 2006;Syed et al., 2009  -Alternative pathways (direct groundwater flows, floodplain inundation, deltaic regions, etc.) are not accounted for.
-Data sharing remains often difficult because of economic and geopolitical constraints and the density of gauge network is decreasing (Shiklomanov et al., 2002).
To counterpart these drawbacks, some authors proposed to use hydrological models (with associated errors and uncertainties) rather than or in combination with in situ data (Trenberth et al., 2007;Dai et al., 2009;Alkama et al., 2011;Haddeland et al., 2011).
An alternative method consists in solving the water budget at the global scale, while distinguishing three main compartments: land, ocean and atmosphere.The latter may be separated into the region over land (land atmosphere) and its complementary part (ocean atmosphere).Figure 1 schematizes water stocks and fluxes involved in the global water cycle.The following equations describe the land Eq. ( 1), ocean Eq. ( 2) and atmosphere Eq. ( 3) global water budgets, respectively (Peixoto and Oort, 1992).
where S l represents the terrestrial water storage, S o the ocean mass, R the global runoff, P the precipitation, E the evapo(transpi)ration, W the total column water vapor and div(Q) the divergence of the vertically integrated water vapor flux.Subscripts l and o designate the spatial average over land and ocean, respectively.The global runoff may be estimated directly from the land or ocean water budgets (Eqs. 1 and 2) (see e.g.Seo et al., 2009;Syed et al., 2010) (ground based and remotely sensed), it may suffer from large uncertainties in some sparsely monitored regions.Moreover, very few direct measurements of evapotranspiration exist (flux towers) and modeling E at the global scale is subject to large uncertainties (Vinukollu et al., 2011).The atmospheric water budget has been introduced by some authors (e.g.Oki et al., 1995;Oki, 1999;Syed et al., 2009) to overcome difficulties in estimating P and E .Dai and Trenberth (2002) showed that the use of the atmospheric water budget improved model based estimates of global runoff.
For the past decade, the ocean mass and the terrestrial water storage can be provided by the Gravity Recovery and Climate Experiment (GRACE) (Tapley et al., 2004;Wahr et al., 2004) space gravimetry mission, as done for instance by Syed et al. (2009) for S l .Before its launch in 2002, no direct measurement of S l and S o were available and since it was difficult to validate estimates from models, variations in water storage were usually neglected (e.g.Oki, 1999;Dai and Trenberth, 2002) then leading to an estimation of the global runoff directly from the net precipitation (P − E ).Nevertheless Syed et al. (2009) showed the importance of taking this term into account by using GRACE estimates of S l .
GRACE products are available since 2002, which may limit the time span of the study to the last decade.Alternatively, Land Surface Models (LSMs) provide monthly estimations of S l with a satisfactory accuracy at basin to global scales, as shown in the various studies comparing GRACE and LSMs (see Ramillien et al., 2008 for a review of studies prior to 2008).The covered period of LSMs depends on the model, but simulations generally run at least over the last two decades.Besides, the ocean mass S o may be derived from satellite radar altimetry observations which provide estimations of the global mean sea level (GMSL) and from in situ hydrographic data.To derive the ocean mass variations, GMSL has to be corrected from the steric component (effect of temperature and salinity), as done by e.g.Syed et al. (2010).
In this study, we use the coupled land-ocean-atmosphere water budgets to estimate the interannual variability of global runoff.Two estimates are computed: one from the land/atmosphere coupling (Eq.4) and the other from the ocean/atmosphere coupling Figures
Section 2 presents the data sets used in this study.In Sect.3, each data set is compared with independent data to ensure its reliability and to give an idea of its uncertainties.The comparison between both global runoff estimates, in terms of interannual variability, is given in Sect. 4. Our global runoff estimate is also compared with global climate indices (ENSO related SOI, AMO).A discussion on the ocean thermal expansion used in the estimation of the ocean mass is provided in Sect. 5.

Data and models used in this study
In this section, we present the data and models used to compute global runoff by the two methods (R l and R o ).Data used for validation purposes are also presented.

Altimetry-based sea level data
For the altimetry-based sea level data, we use the DT-MSLA "Ref" series provided by Collecte Localisation Satellite (CLS; http://www.aviso.oceanobs.com/en/data/products/sea-surface-height-products/global/msla/index.html).This data set is used over the • N to obtain a GMSL time series.The data are further averaged on a monthly basis.

Steric data
Steric sea level is estimated using the updated in situ ocean temperature and salinity data from Ishii and Kimoto (2009), v6.12 (called hereafter IK09).The IK09 temperature data are corrected for the XBT depth bias.The temperature and salinity data are available at monthly interval over 16 depth levels ranging from the ocean surface down to 700 m depth, on a global 1 • × 1 • grid from 1955 to 2009.Steric sea level anomalies are computed over the 0-700 m depth range for the period January 1993 to December 2009.The deep ocean contribution cannot be accounted for since hydrographical data below 700 m are too sparse.Recent studies have shown that in terms of trend the deep ocean contributes by ∼10 % to the total steric effect (Church et al., 2011).Besides, almost all interannual variability in steric sea level is confined in the upper 300-500 m of the ocean (e.g.Llovel et al., 2011).At global scale, salinity does not contribute to the GMSL and is therefore neglected in the following.

Land surface models
To in July 2007.Discrepancies among S l derived from the different LSMs may come from differences in the numerical schemes and meteorological forcing.As noted by Syed et al. (2009), discrepancies among models outputs provide an estimation of the model uncertainties.
Independently, the evapotranspiration modeled by the three LSMs is used for the validation of net precipitation over land (P l − E l ) computed from the atmospheric water budget (Eq.3).

GRACE data
To  ).At each grid mesh, the volume anomalies are obtained by removing the temporal mean.The GRGS data are stabilized during the generation process so that no smoothing or filtering is necessary.When GRACE data is used at basin scales, it has to be corrected from leakage effects due to its low resolution (Longuevergne et al., 2010).Nevertheless such effects have minor impact at the global scale and GRACE data is not corrected here.
As said previously, GRACE provides reliable estimations of spatiotemporal water

Meteorological data
Data used in this study to compute P − E from the atmospheric water budget (W and div(Q)) are based on reanalysis products from the European Centre for Medium-Range Forecast (ECMWF) ERA-Interim data set (Simmons et al., 2007).These are daily global data provided on 1.5 • × 1.5 • grids in units of mm day −1 .All gridded data are further expressed in terms of monthly averages over the period 1993-2009.

Data processing
As the focus of this study is the interannual variability of the global runoff, the seasonal component of each signal presented in the following is removed.This component is obtained by fitting two sinusoidal signals of periods 6 and 12 months.Moreover, although the mean value of the global runoff is still subject to discussions (see e.g.Syed et al., 2009, and references therein), it is not in the scope of this paper.Hence, the temporal mean is also removed from the global runoff estimates.In Eqs. ( 5) and ( 6 As a spherical harmonics (SH) truncation at degree 50 (resolution of 400 km) is applied on GRACE data to obtain water mass variations, we applied the same SH truncation to LSM outputs for a more relevant comparison between GRACE and LSMs derived S l (as suggested by many authors, e.g.Longuevergne et al., 2010).SH truncation is applied only for the model validation, not for the runoff computation.
Finally, all graphs showing temporal evolution of the spatial mean of any variable has been smoothed using a 3-months moving average.

Land/ocean masks: estimate of the high latitudes contribution
As indicated above, altimetry products used here are available only in the 65 • S-65 • N domain.Moreover, ice sheets (Greenland and Antarctica) are generally not modeled in LSMs because of their very specific hydrological behavior, all the more so as very few in situ data are available in these regions and the models validation is then quite difficult.Hence ice sheets and high latitude oceans are excluded from the present study.Figure 2a shows the land and ocean regions considered here.
The exclusion of ice sheets and high latitude oceans has no major consequences in the following since these regions only play a minor role in the interannual variability of global runoff.To assess this, Fig. 2b,c present S l and S o variations derived from GRACE over the four regions shown in Fig. 2a.Despite significant trends in water mass variations -not shown in the graphs -in ice sheets (−160 km 3 yr −1 ) and high latitude oceans (−66 km 3 yr −1 ), these regions are scarcely involved in the global water cycle in terms of interannual variability.the GRACE period has been removed for the purpose of this comparison).The good agreement between the independent LSMs and GRACE derived S l (correlation coefficient of 0.66) reinforces the reliability of both estimates.Figure 3b presents the spatial distribution of the RMS differences between both estimates.Main differences are localized in the inter-tropical zone and more specifically in the Amazon and Congo basins.

Comparison of terrestrial water storage from GRACE and LSMs
Two main reasons may explain this: Eq. ( 1) the hydrological cycle has greater amplitude in these basins than in others leading to higher RMS and Eq. ( 2) GRACE errors are larger near the equator than in high latitudes (Swenson and Wahr, 2006).High RMS between GRACE and LSMs are also found in glaciers regions (Alaska, Scandinavia, Himalaya) which may come from the fact that, contrarily to GRACE, LSMs generally do not account for glaciers.These discrepancies are of less importance in terms of volume since pixel area tends to 0 in high latitudes.
One may note that the depression in S l in 1997-1998 corresponds to a major El-Ni ño event namely characterized by lower/higher precipitation intensity over land/oceans (Gu et al., 2007).

Ocean mass estimate
Global mean sea level (GMSL) variations are the result of variations in the ocean mass (S o ) and in the steric component.Salinity has little effect at the global scale so that the remaining steric component corresponds to the ocean thermal expansion (TE).Ocean mass variations are then derived from GMSL variations corrected from TE according to Eq. ( 6). for the common period [2002][2003][2004][2005][2006][2007][2008][2009].In particular, the TE correction leads to a great negative peak in 2003-2004 that is not shown in GRACE S o .Considering the mitigated efficiency of the TE correction, we prefer not to apply it for the runoff computation done in Sect. 4. Nevertheless, we will come back to this issue in Sect. 5.

P − E over land and atmospheric water budget
Figure 5a shows net precipitation over land (P l − E l ) computed from the atmospheric water budget (Eq.3) and ERA-Interim data.This graph shows a sudden increase in terms of interannual variability (Fig. 6).Besides, Fig. 7 represents the interannual variations of div(Q) and dW /dt over oceans and land.Note that the discontinuity correction has been applied in div(Q) over land and oceans.First, it is shown that dW /dt plays a minor role in the interannual variability of the global water cycle.Second, div(Q l ) is very similar to −div(Q o ) which means that very few water is horizontally exchanged between high-and mid-latitude regions of the atmosphere (separated by the 65 • S and 65 • N parallels).These considerations lead to P l − E l = −(P o − E o ).In that sense, the same offset is applied in the estimation of As said previously and shown in Fig. 7, the 1997-1998 El Ni ño event is characterized by more/less precipitation than usual (not compensated by an increase/decrease in evaporation) over ocean/land.

Global runoff estimates
The comparison between estimates of the global runoff from the land-atmosphere and the ocean-atmosphere water budgets is shown in Fig. 8.Both estimates are very well correlated (correlation coefficient of 0.73, mean standard deviation of 63 km 3 month −1 ).
Interannual variations of global runoff ranges from about −200 to 200 km 3 month −1 , with higher peaks during major ENSO events.In particular, negative peaks in 1994-1995, 1997-1998 and 2009-2010 are related to lower than normal precipitation over land during El Ni ño events, whereas the positive peak in 1999-2000 is related to higher than normal land precipitation during La Ni ña event.
Besides, no significant trend is observed over the whole period.Indeed, the trend equals 0.4 km 3 month −2 for R l and 0.6 km 3 month −2 for R o , which is negligible compared to the mean standard deviation of R l and R o (169 and 228 km  (Dai et al., 2009;Milliman et al., 2008) which found insignificant trends.The difference with Syed et al. ( 2010) may be partly explained by the fact that the authors computed the trend over 1995-2006 and that the global runoff presents higher values before 1995 and lower values after 2006, leading to a reduction of the global trend.Other authors (Dai et al., 2009;Alkama et al., 2011) found non negligible trends over some of the major river basins for the last decades, but these trends seem to compensate against each other at the global scale.
The differences between R l and R o may be explained by two main factors: modeling errors and TE not accounted for in the ocean mass estimation.Concerning the former, considering three different LSMs help us to quantify modeling errors.We find a mean standard deviation of 43 km 3 month −1 due to LSMs discrepancies (represented by the red shading in Fig. 8), which is very low compared to the interannual variability of R.This suggests that the differences between R l and R o are mainly explained by the fact that the ocean mass is not corrected from TE.Nevertheless, the good correlation between R o and R l shows that TE plays a minor role in the interannual variability of global runoff.In Sect.5, we investigate a new method to estimate the ocean thermal expansion component by using the difference R o − R l .

Zonal runoff and correlation with climate indices
As said previously, the global runoff seems to be highly impacted by major ENSO events.Many studies, including Dai et al. (2009), Syed et al. (2010) and Labat (2010), also showed this particular link.In addition, Gu et al. (2007) reported high correlations between ENSO and land precipitation in the intertropical zone and lower but not null correlations in mid-and high-latitudes.In this section, we propose to investigate the correlation between the zonal runoff and different global scale climate phenomena.
First, the zonal runoff is computed over four latitudinal bands separated by the 60 • N, 20 • N and 20 • S parallels (solid lines in Fig. 9).As expected, while the inter-tropical zone contributes for the most part to the interannual variability of the global runoff, northern high latitude and southern mid latitude zonal runoffs are quite negligible.Introduction

Conclusions References
Tables Figures

Back Close
Full For each zone, the zonal runoff is then compared with different linear combinations of several climate indices.The indices considered here are the following: the Multivariate ENSO Index (MEI), the Atlantic Multidecadal Oscillation (AMO), the Pacific Decadal Oscillation (PDO), the Pacific-North America teleconnection (PNA) and the Artic Oscillation index (AO).The reader may refer to Rossi et al. (2011) for a detailed presentation of these indices.For each subset of one or two indices, the linear combination is optimized by maximizing the correlation between the zonal runoff and the combined index.Subsets that give the best results are used to compute a new index called Combined Runoff Index (CRI).Optimization results are given, for each zonal band, by Eq. ( 7 WGHM which is only available since 1992), we also compared CRI and zonal R l over the period 1980-1992, period not used in the calibration of CRI. Figure 9 shows (dashed lines) results of calibration (period 1993-2009) and validation (period 1980-1992); CRI has been normalized to match the range of zonal runoffs variability.The two numbers in the brackets in the legend represent the correlation between zonal R l and CRI for the calibration period and for the overall period, respectively.Figure 9 clearly shows a very good correlation for the inter-tropical zone (namely during the validation period) and more contrasted correlations in mid-and high-latitudes.Not surprisingly, MEI contributes for the most part in the inter-tropical and northern mid-latitude zonal runoffs.For these two zonal runoffs, AMO also plays an important role.Figure 9 and Eq. ( 7) show that these two zonal runoffs are highly anti-correlated, namely with higher/lower runoffs than normal in mid-latitude/inter-tropical zones during Figures

Back Close
Full  (e.g. in 1989 or in 2000).This suggests that during El-Ni ño events, while water mass is shifted westward from the South-American continent to the tropical Pacific Ocean (Gu et al., 2007), it is also shifted northward to mid-latitude continents.Besides, ENSO seems to play a less important role in northern high latitudes and southern mid latitudes.Northern high latitude zonal runoff is logically governed by northern mid to high latitude climate phenomena (Pacific-North America teleconnection and Arctic Oscillation).For each of the three other zones, CRI is a combination of a climate index related to the Pacific Ocean and another related to the Atlantic Ocean.Further investigations are suggested to complete these preliminary results, namely about the relationship between the zonal runoff and climate indices characteristics in the frequency domain (Rossi et al., 2011).

Reconstruction of the ocean thermal expansion
In this subsidiary section, we come back to the aforementioned problem of the ocean thermal expansion (TE) correction and propose an innovative method to reconstruct this component of the global mean sea level.Indeed, as said previously, the observed difference between R l and R o may be mainly attributed to TE.Assuming that Here we propose to use this property to give a new estimation of TE from Eq. ( 8).

Back Close
Full for IK09.For each data base, we computed the thermosteric sea level on a 1 • ×1 • grid at monthly interval since 1993, integrating temperature anomalies from the surface down to 700 m.For that purpose, we first computed density anomalies at each standard level down to 700 m by considering temperature anomalies and using the classical equation of state of the ocean.Then, we integrated density anomalies at each grid point and each time step to obtain the thermosteric sea level (Gill, 1982).
Figure 10 compares TE obtained from our reconstruction and from IK09 (a) and WOD09 (b).
Noting a 3 months lag between our reconstruction and the two other curves -which cause is still under investigation -we have accounted for this delay in Fig. 10.Except at the beginning and end of the time span and around the period 2003-2004, both curves are in quite good agreement in terms of interannual variability, which gives confidence in our method to estimate TE.One may note that the differences between our reconstruction and each dataset are of the same order as the difference between IK09 and WOD09.Concerning the period 2003-2004, this comparison seems to confirm an overestimation of TE derived from IK09 and WOD09 already noted in Sect.3.4 and Fig. 4. To corroborate or invalidate this overestimation, we are presently performing a specific study about the origin of such discrepancies over this period (comparison with TE derived from global circulation models).

Conclusion
The impacts of climate change and anthropogenic factors on the global water cycle represent a critical and timely issue.Namely, the question about an intensification of the global water cycle is highly debated and the answers are still controversial.Evidence of such intensification may be derived by looking at different parameters of the water cycle, e.g. an increase in global runoff (R), with potential implications on basin scale water management.Here we investigated two methods to estimate R over the period 1993-2009.Both methods are derived from the coupling of the land-atmosphere Introduction

Conclusions References
Tables Figures

Back Close
Full and ocean-atmosphere water budgets, respectively.Independent datasets are used to estimate water storage variations in each compartment: atmospheric reanalyses, land surface models, satellite altimetry and in situ ocean temperature data.The two main results are: both estimates correlate very well over the study time span (correlation coefficient of 0.73), giving confidence in the method; no significant trend is observed over the whole period (the trend ranges from 0.4 to 0.6 km 3 month −2 , which is insignificant compared to the mean standard deviation of R ranging from 169 to 228 km 3 month −1 ).
Results also show that major ENSO events largely impact the global runoff (decrease/increase during El Ni ño/La Ni ña).To infer the link between global runoff and climate variability, we computed the zonal runoff over 4 latitudinal bands (intertropical zone, northern and southern mid-latitudes and northern high latitudes) and compared it with different climate indices.We then set up a Combined Runoff Index (CRI) obtained for each zone by optimization of climate indices linear combinations.We showed that CRI well correlates with the zonal runoff even over a validation period not used for the optimization, especially for the intertropical zone which contributes for the most part in the global runoff.Besides, linear coefficients resulting from the optimization provide information about which large-scale climate phenomena are the main drivers for each zonal runoff.Namely, the intertropical and northern mid-latitude zonal runoffs are mainly driven by ENSO and AMO.Moreover, the zonal runoffs in these two zones are anti-correlated, particularly during ENSO events.This suggests that water masses are shifted northward/southward during major El Ni ño/La Ni ña events.Lastly, the difference between the two runoff estimates may be mainly explained by the ocean thermal expansion (TE) term which has not been accounted for in the ocean mass estimate.Consequently, we used our runoff reconstruction to propose an innovative method for providing a new estimate of TE.The comparison with two existing datasets shows a quite good agreement in terms of interannual variability, showing Introduction

Conclusions References
Tables Figures

Back Close
Full the relevance of the method, but point out a period (2003)(2004) during which ocean thermal expansion data may be in error.Introduction

Conclusions References
Tables Figures

Back Close
Full Discussion Paper | Discussion Paper | Discussion Paper | for detailed discussions), among which are the following: -Many regions remain unmonitored.For instance, Fekete et al. (2002) provided global runoff estimates from the 663 major river basins, which represent 71 % of the global runoff.Milliman et al. (2008) used data from 131 river basins, representing 51 % of the global runoff.-The time periods covered by gauge stations are very irregular in terms of start and end time and gaps in observations.Discussion Paper | Discussion Paper | Discussion Paper | , then involving estimations of P and E .Even though precipitation data benefit from dense observation networks Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | estimate the terrestrial water storage component (S l ) we use monthly gridded outputs of three different LSMs: (1) the WaterGAP Global Hydrology Model (WGHM) (D öll et al., 2003); (2) the Interactions between Soil, Biosphere and Atmosphere -Total Runoff Integrating Pathways (ISBA-TRIP) model (Alkama et al., 2011; Decharme et al., 2010); (3) the Land Dynamics (LaD) model (Milly and Shmakin, 2002).WGHM outputs are available through December 2009 whereas ISBA ends in December 2008 and LaD Discussion Paper | Discussion Paper | Discussion Paper | volume variations since 2002 and it would have been possible to complete GRACE data with LSMs outputs and altimetry based ocean mass for the period 1993-2002.In order to keep consistency in our computations over 1993-2009, we prefer not to use GRACE products to estimate terrestrial water storage and ocean mass variations.Nevertheless, in both cases, comparison with GRACE observations is performed to increase confidence in data.Discussion Paper | Discussion Paper | Discussion Paper | ), S l and S o are derivated with respect to time and any trend in S l and S o would lead to constants in the runoff.Consequently, S l and S o are detrended in the following.Discussion Paper | Discussion Paper | Discussion Paper |

Figure
Figure 3a shows the total water storage derived from the LSMs; the blue shading represents the root mean square (RMS) deviation of each model with respect to the average.The mean standard deviation (RMS) between LSMs is 1.26 mm which is quite low compared to the amplitude of the interannual variations of S l .The comparison with

Figure 4
Figure 4 presents the ocean mass variations S o derived from GRACE and computed from GMSL, corrected or not from the thermal expansion TE.The thermal expansion correction over 2002-2005 clearly deteriorates the correlation with GRACE derived S o (correlation coefficient of 0.17 and 0.45, respectively with and without correction, 4643 Figure5ashows net precipitation over land (P l − E l ) computed from the atmospheric water budget (Eq.3) and ERA-Interim data.This graph shows a sudden increase inP l − E l in2003.The difference between the temporal mean over 1993-2002 and the one over 2003-2009 equals 3 mm month −1 , which is of the same order of magnitude as the interannual variability.We are currently investigating the causes of this important shift, but we suspect a change in input data in the ERA-Interim procedure.Namely radiance measurements from AIRS are assimilated since July 2003, which may have an impact on the global water cycle (P.Poli, personal communication).To assess the artificial origin of the discontinuity, we compare P l − E l from ERA-Interim with independent precipitation and evapotranspiration datasets (see Sect. 2).We use the outputs of the three LSMs presented previously (ISBA, WGHM and LaD) to estimate evapotranspiration.The blue curve and blue shading in Fig. 5b show the mean and standard deviation of the 18 computed time series of P l − E l .The observed P l minus modeled E l is very similar to the ERA-Interim estimate before 2003 but it does not present any shift since 2003, which reinforces the assumption of an artificial origin of the discontinuity.We then decided to correct ERA-Interim P l − E l by adding a constant offset over the period 2003-2009.The value of the offset is obtained by minimizing the difference between ERA-Interim P l − E l and observed P l minus modeled E l .The optimum offset value (−2.78 mm month −1 ) leads to a very good correlation between the two estimates.The same analysis may have been done over oceans but very few evaporation data over ocean exist and comparison between two of the main existing data sets (OA-Flux,Yu and Weller, 2007; HOAPS, Andersson et al., 2007)  shows large discrepancies in 3 month −1 , respectively).These values are quite lower than the one found by Syed et al. (2010) (540 km 3 yr −2 or 3.75 km 3 month −2 ) but are in agreement with results of some other previous studies Introduction Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | ) (climate indices have been normalized before the optimization).CRI(−60/ −20) = 0.53 × PDO + 0.47 × AMO CRI(−20/ +20) = −0.62 × MEI − 0.38 × AMO (7) CRI(+20/ +60) = 0.54 × MEI + 0.46 × AMO CRI(+60/ +90) = 0.73 × PNA + 0.27 × AO Since LSMs and ERA-Interim outputs are available since 1980 or earlier (apart from Discussion Paper | Discussion Paper | Discussion Paper | major El-Ni ño events (e.g. in 1983 or in 1998).The reciprocal is true for La-Ni ña events Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Tapley, B. D., Bettadpur, S., Watkins, M. M., and Reigber, C.: The gravity recovery and climate experiment; mission overview and early results, Geophys.Res.Lett., 31, L09607, doi:10.1029/2004GL019920,2004.Trenberth, K. E., Smith, L., Qian, T. T., Dai, A., and Fasullo, J.: Estimates of the global water budget and its annual cycle using observational and model data, J. Hydrometeorol., 8Discussion Paper | Discussion Paper | Discussion Paper | igure 1: Stocks and fluxes in the global water cycle.Refer to the text for a detailed escription of the variables.S l represents the terrestrial water storage, S o the ocean mass, R e global runoff, P the precipitation, E the evapo(transpi)ration, W the total column water apor and div(Q) the divergence of the vertically integrated water vapor flux.Subscripts l and designate the spatial average over land and ocean, respectively.

Fig. 1 .
Fig. 1.Stocks and fluxes in the global water cycle.Refer to the text for a detailed description of the variables.S l represents the terrestrial water storage, S o the ocean mass, R the global runoff, P the precipitation, E the evapo(transpi)ration, W the total column water vapor and div(Q) the divergence of the vertically integrated water vapor flux.Subscripts l and o designate the spatial average over land and ocean, respectively.

Fig. 2 .
Fig. 2. (a) Land and ocean masks (high latitude Arctic and Antarctic oceans in light blue and ice sheets in green).(b) Detrended interannual variations of GRACE derived water stocks in land (except ice sheets) and oceans (except high latitudes).(c) Detrended interannual variations of GRACE derived water stocks in ice sheets and high latitude oceans.

Fig. 4 .
Fig. 4. Comparison of global mean sea level (GMSL) computed from altimetry, GMSL corrected from the thermal expansion TE and GRACE ocean mass.

Fig. 5 .
Fig. 5. (a) Interannual variability of precipitation minus evapotranspiration over land (P l − E l ) computed from the atmospheric water budget and ERA-Interim data.(b) Corrected P l − E l with a uniform offset of −2.78 mm month −1 over 2003-2009, and observed P l minus modeled E l used to find the optimum offset (blue shading represents discrepancies between P l −E l obtained from the different the data sets).