Radar Altimetry as a Proxy for Determining Terrestrial Water Storage Variability in Tropical Basins

Abstract

The Gravity Recovery and Climate Experiment (GRACE) mission has provided us with unforeseen information on terrestrial water-storage (TWS) variability, contributing to our understanding of global hydrological processes, including hydrological extreme events and anthropogenic impacts on water storage. Attempts to decompose GRACE-based TWS signals into its different water storage layers, i.e., surface water storage (SWS), soil moisture, groundwater and snow, have shown that SWS is a principal component, particularly in the tropics, where major rivers flow over arid regions at high latitudes. Here, we demonstrate that water levels, measured with radar altimeters at a limited number of locations, can be used to reconstruct gridded GRACE-based TWS signals in the Amazon basin, at spatial resolutions ranging from 0.5 to 3°, with mean absolute errors (MAE) as low as 2.5 cm and correlations as high as 0.98. We show that, at 3° spatial resolution, spatially-distributed TWS time series can be precisely reconstructed with as few as 41 water-level time series located within the basin. The proposed approach is competitive when compared to existing TWS estimates derived from physically based and computationally expensive methods. Also, a validation experiment indicates that TWS estimates can be extrapolated to periods beyond that of the model regression with low errors. The approach is robust, based on regression models and interpolation techniques, and offers a new possibility to reproduce spatially and temporally distributed TWS that could be used to fill inter-mission gaps and to extend GRACE-based TWS time series beyond its timespan.

1. Introduction

The Gravity Recovery and Climate Experiment (GRACE) satellites [1] have provided unprecedented information on global water-storage dynamics by estimating variations in the terrestrial water-storage (TWS) variability. GRACE data has been used to evaluate global TWS trends [2], to monitor and characterize droughts (e.g., [3,4,5], estimate groundwater depletion due to anthropogenic activities [6,7,8,9] and improve hydrological models through data assimilation techniques (e.g., [10,11]) towards a better quantification of land water-storage components and their spatial and temporal variability. Since GRACE data provides estimates of the variability of the total amount of water-storage change within a domain, the inability to identify the water layer and state has been a challenge in decomposing TWS signals into its different components (i.e., groundwater, soil moisture, snow and surface water (including rivers, floodplains, wetlands, lakes and reservoirs). Few studies attempted to determine the impact of surface water storage (SWS) on TWS at the basin scale (e.g., [5], [6]). Using state-of-the-art hydrological models, [12] were able to decompose TWS and determine the spatially-distributed impact of its components. Their results indicate that, although soil moisture has the major impact on global TWS, surface water has a key role in the tropics and over large rivers. In the Amazon basin, for example, SWS accounts for as much as 27% of TWS variability. Other studies, based on remote-sensing data using different metrics, revealed that SWS amplitude is about 50–60% of that of TWS (e.g., [13,14,15]). As well, using in situ data, [16] found a high correlation between river dynamics and TWS in the Amazon basin.

GRACE data provides a 15-year span, from April 2002 to June 2017, with gaps. Its relatively short temporal availability limits its applicability to long-term analyses, but the recent launch of the GRACE Follow-On (GRACE-FO) mission, in May 2018, will prolong global TWS time series for several more years. In order to fill the one-year gap between missions, and to extend TWS time series beyond its timespan, different techniques estimating TWS variability based on model outputs [17], [18], in situ data [19], the combination of both [20], and reanalyses [21] have been presented in the literature.

Table 1. Previous attempts to reconstruct or predict terrestrial water-storage (TWS) signal.

Reference	Approach	Domain	Time span	Gridded
[17]	Water budgetbased on model output	Basins in Europe and USA	2002–2006	No
[19]	Water budget based on observed data and model output	Amazon and Mississippi basins	1970–2006 (Amazonas) 1928–2006 (Mississippi)	No
[41]	Water budget based on observed and reanalysis data	Basins in several continents	1992–2005	No
[42]	Reduced optimal interpolation using GRACE data and in situ river level records	Amazon basin	1980–2008	Yes
[18]	Water budget using GRACE data and model outputs	Amazon basin	1948–2010	No
[20]	Statistical model using reanalysis data	Global	1985–2015	Yes
[43]	Physically based modelling and deep learning	India	2002–2017	Yes

provides a list of these studies. Although TWS derived from these techniques showed an overall agreement with GRACE data, limitations in representing TWS variability and trends by models [22,23] and restricted access to ground-based observations could hinder their widespread application. In this sense, an alternative solution to fill these gaps is the use of other satellite sensors as a proxy to determine TWS.

Radar altimetry (RA) technology has been continuously used to measure water levels in oceans and inland waters since 1992 (e.g., [24,25,26]). Numerous studies have shown its potential (i) to monitor large-scale inland water levels (e.g., [26,27,28,29,30], floodplain dynamics and water-storage change [13,31,32], to estimate river flow velocity and streamflow [15,33,34], and (iii) to calibrate [35,36,37] and evaluate hydrological models [12,38]. Here, we explore the relationship between GRACE data and surface water-level dynamics, and the potential to derive TWS from large RA datasets. Based on previously describe findings, our hypothesis is that, over tropical basins, where surface water significantly contributes to TWS variability, radar altimetry can be used to reproduce terrestrial water-storage signals. In that regard, we used a RA network composed of 984 virtual stations (VS; i.e., the intersection between satellite tracks and water bodies) in the Amazon basin and developed an approach based on regression models and interpolation techniques to transform point-based RA signals into a gridded TWS product. In order to fill out areas without RA data availability, we evaluate three interpolation techniques (kriging, inverse distance weighting and nearest neighbor) and, based on a sensitivity analysis, we estimate the errors related to each technique and the maximum search radius for interpolation. Calibration and validation periods are defined in order to determine the potential for the technique to be extrapolated to other periods where GRACE-based TWS estimates are not available. At the basin scale, surface water, soil moisture and groundwater in tropical regions are climatologically in phase, as shown in [12]. This can be explained by the interconnectivity between surface and soil waters, as suggested in previous studies [39,40]. In this sense, by exploring the relationship between surface water and TWS, we assume that the SWS signal is a major contributor to TWS variability in the Amazon basin and that its seasonality is generally correlated to the other hydrological variables.

2. Materials and Methods

2.1. Gravity Recovery and Climate Experiment (GRACE) Terrestrial Water Storage (TWS)

GRACE measurements of the Earth’s time-variable gravity are processed by different centers that provide monthly gravity solutions using two approaches: spherical harmonics and mass concentration (mascons) blocks [44,45]. During the first years of GRACE observations, spherical harmonics was the standard. However, recent studies have shown that the mascons solution offers several key advantages, including (1) reduced leakage through the land-ocean interface, which increases signal amplitude; (2) removal of north-south stripes is achieved without the need of empirical filters; and (3) finer resolution for smaller areas [45,46]. Here, we used the University of Texas’ Center for Space Research (CSR) version RL05 [47]. Although CSR GRACE masons are natively estimated at 120 km-wide mascon blocks and resampled to a 0.5° spatial resolution grid, the actual resolution of its solution is about 250–300 km along the equator line. Mascons solutions from other GRACE data-processing centers have different native resolution (e.g., three degree equal-area caps). Therefore, our analysis considered four grid sizes: 0.5, 1, 2 and 3 degrees. For this reason, TWS was resampled into coarser resolutions by averaging the 0.5° data. GRACE data have multiple sources of errors that emanate from instruments measuring through post-processing [48]. A number of attempts to assess GRACE estimates have been undertaken (e.g., [49,50]). GRACE data accuracy has been estimated to be less than 1cm in equivalent water height, when averaged over areas larger than about 4 × 105 km², and errors increase as the area under observation decreases [51]. [49] noted that the errors increase at lower latitudes (~2.6 cm), whereas it is around 0.8 cm near the poles. While evaluating the global CSR RL05 mascons solution, [46] noted that GRACE-based TWS errors were the lowest (<3cm) for large basins (>0.5 million km²), but up to three times higher for small basins.

2.2. Radar Altimetry Data

A total of 984 virtual stations are available over the Amazon basin for 2002–2016 on the Theia-Land’s Hydroweb platform [29], derived from multiple satellite missions: Jason (11 VS), Jason 2 (88 VS), Envisat (741 VS) and Sentinel 3A (144). The temporal resolution varies as a function of the satellite, i.e., 10 days for Jason, 27 days for Sentinel 3A and 35 days for Envisat. [28] list two main sources of errors in RA data. The first can be due to bounced back radar echoes, captured by the antenna receiver that are caused by solid reflectors (e.g., vegetation, islets, etc). When such echoes are combined with those from water surfaces at a different level, the height inference will be affected. Second, RA retrievals can be biased if the detected surface is off-nadir, because the signal received by the satellite may be dominated by water surface at the edge of the footprint. A tool was developed by [28] to overcome those limitations and identified errors in altimetric time series below 20 cm up to meters. More specifically, they reported a characteristic error in the Envisat series below 30 cm. More details about the generation of altimetry-based water levels used in this study can be found in [28].

2.3. Estimating Terrestrial Water Storage from Radar Altimetry

In this study, water level data was filtered by removing VS containing less than 12 months of data. Gaps in the VS time series were filled using a shape-preserving piecewise cubic interpolation [52]. This method requires at least four neighboring grid points for the interpolation to be performed. To reduce the chances of generating unrealistic interpolations, we only interpolated time series whose gaps do not exceed 30% of its length. For example, a time series spanning 24 months was gap filled only if the number of valid records exceeded 17, i.e., a maximum of 7 gaps (Figure S1 shows an example of RA time-series interpolation using this technique). Gaps in the tails of the time series were not filled, i.e., we did not perform temporal extrapolation. The filtering process resulted in 839 VS, whose availability is unevenly distributed through time, i.e., a fraction of that figure is available each month (Figure S2 shows the temporal distribution of VS data availability). For example, from 2004 to 2008, we used ~600 VS; from 2008 to 2010, ~700 VS were available. Because few stations are available after 2010, we limited the study period from October 2002 to July 2010, which was divided in two parts: five years for the calibration period (April 2002–March 2007) and the remaining three years for the validation period (April 2007–July 2010).

Radar altimetry data was converted to the monthly time step to match GRACE-based TWS estimates (TWS_GC) by averaging observations within the same month. GRACE pixels corresponding to each VS were identified and correlations were computed between water levels and TWS. Because some RA-based water level anomalies correlate better with TWS than others, we selected virtual stations with highest correlation (Figure 1). This procedure was repeated for each of the four GRACE resolutions considered in this study. We then used the best correlated points to fit a linear (TWS = a·H + b, where H stands for water level anomalies) and a quadratic (TWS = a·H² + b·H + c) regression models that estimate TWS from radar altimetry (Table S1 provides regression equations and corresponding metrics for the 3° solution). For each GRACE pixel, the regression model with the lowest mean square error was selected. The regression models obtained for the calibration period were used to derive the estimates for the validation period, i.e., during validation, the procedures described above are not repeated. As a result, the number of pixels with regression-based TWS (TWS_reg) estimates varied according to the spatial resolution and month (i.e., finer resolutions result in more pixels covering the domain, hence more VS are used, one for each pixel, wherever data is available) (Figure S3). The average number of GRACE pixels with water-level data for each resolution was: 388 (0.5°), 216 (1°), 81 (2°) and 41 (3°).

Note that TWS_reg estimates are limited to the pixels containing at least one VS. The remaining pixels do not have RA data with which a correlation could be established. To generate interpolation-based TWS estimates (TWS_interp) in those remaining pixels, three interpolation methods were tested: ordinary kriging (krig), inverse distance weighting (IDW), and nearest neighbor (NN). The method with best performance was applied to produce the final TWS fields. Kriging is a geostatistical method based on the regionalized variable theory, which accounts for the spatial variability of natural properties by means of a stochastic approach [53]. It has been widely used to perform spatial analysis of meteorological and hydrological data (e.g., [54,55,56,57]). Most of these data often present a spatial dependence, i.e., similarities between two points are likely to increase as the distance between them decrease. Thus, one can infer how a given property varies in space between sample points based on an experimental variogram, which describes the spatial variability of a regionalized property or variable.

Here, we used the time-averaged TWS from GRACE to derive a unique experimental variogram for each spatial resolution, from which weighting factors (λ_i), associated with TWS_reg at location x_i, are derived. The kriging estimate TWS_interp at location x₀ is given by:

$$TW{S}_{interp}(x_o)={\displaystyle {\sum }_{i=1}^n{\lambda }_{i}TW{S}_{reg}}(x_i)$$

(1)

In the IDW method, TWS_interp is the weighted average of the known neighbor values; the weights (w_i) are computed as:

$$w_i=\frac{1}{d_i^2}$$

(2)

where d is the distance between the x_i and x₀. Both krig and IDW methods were applied using the R package gstat [58].

The NN method uses the Voronoi tessellation of a limited number of sample points to select, for each grid point (x_o), neighboring data points for interpolation [59]. Weighting is based on the influence area of such data points. We used Matlab’s internal ‘scatteredInterpolant class’ to perform the NN interpolation. Previous studies have applied IDW and NN to (spatially) interpolate rainfall and other environmental data at multiple time scales (e.g., [60,61,62,63,64,65]).

Finally, the different TWS_interp estimates were merged with TWS_reg to produce the final altimetry-based TWS (TWS_alt). The interpolation method with best performance was applied to produce the final TWS fields.

2.4. Evaluation Procedure

For both calibration and validation periods, the mean absolute error (MAE), correlation coefficient (CC) and standard deviation ratio (γ) between our TWS estimates and TWS_GC were calculated.

$$MAE={\displaystyle {\sum }_{i=1}^n}\frac{\left|TW{S}_{GC}-TW{S}_{est}\right|}{n}$$

(3)

$$\gamma \_=\frac{{\sigma }_{alt}}{{\sigma }_{GC}}$$

(4)

where TWS_est can be TWS_reg, TWS_interp or TWS_alt; σ_GC and σ_est are, respectively, standard deviation of TWS_GC and TWS_est time series.

We also evaluated the final product during extreme events and by comparison against other products. Regarding the former, we analyzed TWS_GC seasonal anomalies between 2003 and 2010; and selected two extreme periods (one wet and one dry) to carry out such comparison. We then compared TWS estimates averaged for the entire Amazon basin and its main tributaries: Solimões, Negro, Tajapós, Madeira, Xingu rivers, and the area where the extreme events occurred.

Finally, we compared our basin-averaged estimates against previous studies that attempted to reconstruct TWS. In particular, we used the estimates produced by [12] as a benchmark for our product. Such TWS estimates are a result of state-of-the-art hydrological models, i.e., the Noah land surface model with multi-physics options (Noah-MP: [66] and the Hydrological Modeling and Analysis Platform (HyMAP) river-routing scheme [67,68], simulating five major water storage components (i.e., groundwater, soil moisture, snow, canopy interception and surface water storage, which includes rivers and floodplains) at the global scale. Anthropogenic activities were neglected in their TWS estimate. The final product is the mean of a 12-member ensemble composed of model runs forced with different meteorological and precipitation datasets. More details on the ensemble generation can be found in [12].

3. Results

3.1. GRACE-Based TWS versus Radar Altimetry

Radar altimetry data correlates well with TWS_GC, with CC⩾0.8 in most pixels (1st row in Figure 2). CC tends to increase as the spatial resolution degrades. At finer resolutions (0.5° and 1°), high correlations (CC > 0.8) were mostly found in the central Amazon basin. Few points exhibited a weak, although positive, relationship (CC < 0.4) between TWS_GC and RA data, and are located near the basin’s boundary.

3.2. Evaluation of the Regression-Based TWS Estimates

Supported by the good correlation between TWS_GC and water-level observations in the study area, we estimated TWS from RA data. TWS_alt and TWS_GC correlates well and deviations measured by the MAE indicate that errors are low (2^nd and 3^rd rows in Figure 2). Previous studies (e.g., [19]) showed Amazon basin-averaged TWS amplitudes of 15 cm, with values above 170cm over the center of the basin (see Figure S4). Few pixels showed MAE >10cm, and even less reached the maximum error (MAE ~20cm). Most estimates had MAE lower than 10 cm, corresponding to <10% of local TWS_GC amplitude.

Following the relationship noted between GRACE-based TWS and RA data, our estimates were closer to TWS at coarser resolutions. At 0.5° and 1°, the averaged CC were similar (~0.78), although MAE values were ~0.4 cm lower at 1° than those obtained at 0.5°. At 2 and 3°, CC values were mostly higher than 0.8 and MAE values lower than 5.6cm. Considering the straightforwardness of the proposed approach, the errors in most pixels were relatively low, indicating that the estimates were satisfactory. Evaluation metrics did not show much change between calibration and validation periods (

Table 2. Amazon basin-averaged performance of regression-based TWS (MAE = Mean Absolute Error; CC = Correlation Coefficient; γ = standard deviation ratio).

Product—Resolution	MAE (cm) (calib./valid)	CC (calib./valid)	γ (calib./valid)
TWSalt−0.5°	6.43/6.35	0.78/0.79	0.77/0.77
TWSalt−1°	5.9/6.05	0.78/0.78	0.78/0.75
TWSalt−2°	5.41/5.58	0.82/0.82	0.81/0.79
TWSalt−3°	5.09/5.40	0.80/0.78	0.79/0.75

).

3.3. Interpolation Evaluation

TWS_interp was satisfactory at all spatial resolutions. Comparisons against TWS_GC show high CC (>0.8) and low MAE (<10cm) at most interpolated pixels (Figure 3). Correlation was lower for pixels near the basin boundary because of the lower number of VS available within a short radius, reducing the quality of the interpolation. TWS_interp also shows increased MAE values towards the basin outlet, especially around the mainstream in the eastern part of region, where TWS_alt also exhibited high MAE (>15 cm). Given the high RA data availability in that part of the basin, the low performance of both TWS_alt and TWS_interp could be related to the sea tide influence on the water level. According to [69], the tide effect can be noticed more than 1000 km upstream from the Amazon basin outlet. This means that the monthly mean water level derived from altimetry could be contaminated by the tide, reducing the correspondence with GRACE-based TWS.

The interpolation results in Figure 3 were obtained by setting the maximum search radius (MSR) varying from 3.5° to 4.5°, depending on the grid resolution, and a minimum of three pixels with VS within that radius. To test how the search radius affects the interpolation, we calculated the mean CC resulting from the kriging and NN interpolation as such radius increases. This analysis is necessary because the quality of the estimates decreases as MSR increases, but setting a small MSR means neglecting pixels near the basin boundary. Because the results from IDW and kriging methods were similar, Figure 4 shows the plots for kriging and nearest neighbor (NN) methods only. We noted that, using kriging (and IDW), correlation is above 0.8 for all MSR and about 100% of pixels were filled with MSR ≥4.5° for 3° spatial resolution and MSR ≥3° for finer spatial resolutions. For the NN method, the gain of spatial coverage comes with a high cost to the quality: CC decreases faster than that using kriging as MSR increases, with CC falling below 0.8 for MSR ≥2.5° and 10% of empty pixels (Figure 3). Based on that evaluation, we limited the search radius at 3.5° for the 0.5° grid, 4° for the 1° grid, and 4.5° for 2° and 3° grids.

3.4. Basin-Average TWS

We aggregated all altimetry-based TWS estimates (TWS_reg + TWS_interp) and averaged the final product (TWS_alt) over the Amazon River basin and basins of its main tributaries. Among the droughts and floods that have occurred in the region, we selected two extreme events, one during the dry season (MAM 2004) and one wet (DJF 2009), that covered practically the same area extent in the northern part of the Amazon Basin (Figure 5). The time series for such ‘extreme’ area are shown in Figure 5, as well as selected metrics (MAE, CC and γ) between TWS_GC and TWS_alt for each basin and resolution.

TWS_alt compared well with TWS_GC over the study period (October 2002–July 2010), especially at 1° spatial resolution, which was used to plot the time series in Figure 5. The amplitude, phase and shape of the estimated time series are very close to GRACE data. TWS_alt for the ‘Extreme’ area is also similar to TWS_GC, indicating that our estimates are not significantly affected by such extremes and can provide reliable TWS estimates under different scenarios (regular, droughts and floods).

Although the difference between the seasonal anomaly detected by TWS_GC and TWS_alt were relatively large for some pixels during MAM 2004, such deviations did not affect the basin average time series. For the wet period (DJF 2009), the seasonal anomalies are well represented both at the pixel scale as well as in the basin-average time series. As the spatial resolution degrades, we can see three different trends in MAE: an increase in the Negro and Solimões River basins and in the selected extreme event region; a decrease in the Madeira, Tapajós and Xingu River basins; and little change in the Amazon River basin (Figure 6).

Except for Xingu at 0.5 and 1°, and Negro at 3°, MAE between TWS_alt and TWS_GC time series was below 6 cm. With a few exceptions, CC was above 0.95 for all resolutions and basins. The lowest CC was found for the Negro River basin at 3°. At 3° and compared to other basins, TWS_alt for the Negro River basin was not as well represented: MAE >6cm and CC <0.9. This limitation is linked to the fact that, at such resolution, two among the few pixels in the Negro River basin have low CC between RA and TWS_GC (see Figure 1). For finer resolutions, the number of well correlated virtual stations is higher, enhancing the quality of our estimates in that basin. The standard deviation ratio (γ) was the only metric with a similar (negative) trend among the analyzed cases, i.e., it decreased for all sub-basins as the spatial resolution degrades. Regardless of the basins and spatial resolution, γ > 1 show that the amplitude of TWS_alt signal is greater than that from TWS_GC. In most basins, γ decreases as spatial resolution degrades from 1° to 2°. At 3°, γ converged to a range between 1.5 and 1.2 for the Amazon, Madeira, Negro and Tapajós basins, and below 1.1 for Extreme, Solimões and Xingu River basins.

4. Discussion

TWS_alt is competitive compared to previous studies that have reconstructed or predicted TWS time series (e.g., [17,18,21]). As an attempt to reconstruct gridded TWS using the Global Land Data Assimilation System (GLDAS; [70]), [18] found that TWS estimates are poorly correlated with GRACE data, with CC <0.5 in about 50% of the basin. They also proposed a linear regression model to predict the basin-averaged TWS signals that satisfactorily captured the general variation patterns found in TWS_GC (no metric was reported for that product). Based on a visual comparison between their plots and ours confirm the competitivity of our product. Not only our time series depicted the variation patterns of TWS_GC, we show a good agreement for several sub-basins of the Amazon basin in terms of both time series and spatial distribution.

Water budget-based TWS time series by [20] was also derived for the entire Amazon basin. Performance metrics relative to less than a four-year period during which GRACE data was available are not reported in his paper. The largest difference reported between TWS_GC and its estimates occurred in 2003, ranging from five to ten cm. Using seven time series, as shown in Figure 6, each ~eight-year long, we quantified the deviations to the basin-averaged GRACE-based TWS. Only seven months (3%) had deviations larger than 15 cm, whereas 300 months (~47%) had deviations smaller than 5 cm (Figure S5).

Reference [41] found correlations ranging from –0.31 to 0.77 between their reanalysis-based TWS and GRACE estimates for a number of basins across North America, Europe, Asia and Australia. Using water level data at 58 in situ stations, [19] reconstructed TWS signal over the Amazon using in situ water-level data and obtained correlations of about 0.9 between their time series and that from the basin-averaged GRACE data. [21] reconstructed global 0.5° TWS time series over ~30 years and obtained correlations (0.94–0.97) and MAE (2.3 cm–5.1cm) similar to ours, over the Amazon basin and sub-basins. In terms of standard deviation ratio, their product overperformed ours: while we found 1.02 < γ <1.17, that metric ranged from 1.02–1.06 with their dataset.

5. Conclusions

Here, we presented a simple and efficient method to estimate monthly gridded TWS over the Amazon Basin, at different spatial resolutions, fully based on remote sensing data. This study offers a straightforward option for gap-filling GRACE-based TWS data and reconstruct past TWS over regions where surface water accounts for an important part of TWS signal, such as the Amazon basin [6]. Unlike previous approaches that relied on ancillary information (e.g., hydrological model outputs, meteorological forcing and reanalyzes), the proposed method only requires radar altimetry time series. It consists of using regression coefficients between GRACE-based TWS and RA-based water levels, to predict TWS solely from radar altimetry data. Radar altimetry data is available from 1992 to present and can be acquired near real time, allowing a low-latency TWS estimate, as opposed to GRACE’s several-month latency.

We showed that altimetry-based TWS compares well with GRACE-based TWS, as indicated by the performance metrics, and that selected extreme events can be satisfactorily represented. Although altimetric data is not available with enough density to provide TWS estimates all over the domain at finer resolutions (i.e., 0.5° and 1°), we were able to reproduce spatially continuous gridded products, at such grid sizes, by means of interpolation of regression-based TWS estimates. Results also demonstrated that TWS time series can be reconstructed with high accuracy at coarser spatial resolutions over the Amazon basin with as few as 41 virtual stations.

Despite the promising results, future analyses are needed in order to determine the applicability of the method over other regions and, ultimately, at the global scale. The main limitations that need further investigation could be found in regions with low RA data availability or low influence of surface water on TWS, in particular snow-dominated regions. Although the radar altimetry dataset adopted in this study provides a large number of virtual stations over the Amazon basin, data gaps limited its use to the 2002–2010 period. As a consequence, the recent development of tools to automatically extract multi-mission radar altimetry data (e.g., [43]) could contribute to increase data availability. Nonetheless, our approach has a great potential to provide high-temporal and -spatial estimates of TWS in a near future, with the launch of the Surface Water and Ocean Topography (SWOT) mission, scheduled for 2021.