Background & Summary

Runoff and evapotranspiration (ET) are principal elements of the water cycle, markedly influencing the interchange of water, energy, and carbon among the terrestrial biosphere, hydrosphere, and atmosphere1,2,3. Within the global hydrological cycle, surface precipitation is partitioned into runoff (approximately one-third) and ET (about two-thirds). Notably, the ET consumes more than half of the total global absorbed solar radiation4,5,6,7. The Tibetan Plateau (TP), often referred to as the “Asian Water Tower”, is the source of numerous large rivers. It plays a pivotal role in supplying water to both ecological environments and human societies in downstream areas8,9,10. In the context of global change, the TP serves as a potential “tipping point” in the Earth system, exerting substantial regional and global influences11,12,13. Currently, the TP is currently experiencing significant cryosphere ablation owing to regionally amplified warming, which has a marked impact on water, energy, and carbon cycles, and lead to increased occurrences of extreme hydrological events and uncertainties in the carbon budget14,15,16,17,18,19,20,21.

To date, generating reliable runoff and ET datasets for high-altitude river basins across the TP has been challenging, hindered by the scarcity of quality-controlled observational data22,23,24,25,26,27,28. This lack is exacerbated by harsh environmental conditions, infrastructural limitations, and restrictive data-sharing policies, particularly for transboundary river data29. Furthermore, the limited availability of observed data present challenges in accurately depicting the large-scale spatial heterogeneity of hydrometeorology in the high-mountain basins.

Addressing these challenges necessitates the creation of large-scale gridded runoff and ET products across the TP by integrating sparse ground observations with extensive remote sensing and reanalysis datasets. These could achieve through physics-based hydrological models or data-driven machine learning and multi-source data fusion techniques. The water input in the TP’ water balance is supplied not only by precipitation but also by the storage and regulation from glaciers, snow, and permafrost30,31,32; however, numerous current research approaches of gridded runoff and ET products do not adequately account for cryosphere-hydrological processes25,26,27,28,33,34,35,36,37,38,39,40,41,42,43.

Given the TP’s significant elevation gradients and diverse topography, the relatively coarse spatial resolution of existing products fails to accurately capture fine variations in runoff and ET, particularly within montane river valleys11. High-resolution hydrological modeling shows significant advantages in detailed exploration of complex topographic areas, thereby markedly improving the accuracy of hydrological simulations compared to those employing lower resolutions22,44,45,46,47.

This study aims to generate a high-quality 5-km spatial resolution monthly gridded dataset of runoff and ET for the headwater regions of major rivers across the TP, referred to as TPRED. These rivers include the Yellow, Yangtze, Mekong, Salween, Brahmaputra, Ganges, and Indus rivers, covering the TP from east to west, spanning 1998 to 2017. The TPRED was generated using the advanced cryosphere-hydrology WEB-DHM model, which was constrained by valuable observed discharge at the outlet of each headwater. This dataset is positioned to significantly contribute to various research fields including hydrology, meteorology, water-carbon cycle studies, water resource management, and the improvement of hydrological model accuracy, thereby advancing the goals of sustainable development for water resources across the TP.

Methods

Study area

The Tibetan Plateau (TP) constitutes the third-largest cryospheric area globally, following only Antarctica and the Arctic48,49. At an average altitude of over 4,500 meters, the TP is recognized as the highest plateau on Earth. The TP receives its water supply from both monsoon and westerly, alongside the meltwater from glaciers, snow and frozen soil, making it a crucial water source for numerous major rivers that extend from High Mountain Asia to the densely populated lowlands8,9,10,43.

This study focuses on the headwaters of seven major exorheic rivers across the TP from east to west: the Yellow, Yangtze, Mekong, Salween, Brahmaputra, Ganges, and Indus rivers (Fig. 1a). These headwaters exhibit diverse hydrometeorological characteristics, shaped by their individual responses to climatic and environmental changes50,51. Details regarding these river headwaters are presented in Table 1.

Fig. 1
figure 1

Headwater regions of seven major rivers across the TP. (a) The hydrological observation stations located at the outlet of each headwater region are Tangnaihai (Yellow River), Zhimenda (Yangtze River), Changdu (Mekong River), Jiayuqiao (Salween River), Nuxia (Brahmaputra River), Karnali (Ganges River), and Besham (Indus River). (b, c) Annual average precipitation, temperature and landcover types and proportion across seven headwaters. See Table 2 for the sources of these data.

Full size image
Table 1 Basic characteristics of seven headwater regions across the TP.
Full size table

High elevations lead to almost half of the region experiencing average annual temperature below freezing, notably in the headwaters of the Indus and Yangtze (Fig. 1b and Table 2). These areas have shown a sharp increase in temperatures over recent decades, with warming rates significantly exceeding those of comparable latitudes14. Precipitation distribution varies significantly across these headwaters, with higher rainfall observed in the southern regions of the Ganges and the southeastern areas of the Brahmaputra, often exceeding 1,000 mm annually (Fig. 1c and Table 2). Climate change has further altered precipitation patterns, increasing their irregularity and complexity15.

Table 2 Datasets used in model input and evaluation.
Full size table

The distinct climatic and environmental conditions across the TP support a diverse ecosystem with complex vegetation patterns. Grasslands dominate half of the landscape across these headwaters (Fig. 1d and Table 2). A warming and moistening of the climate has induced a visible “greening” effect on the vegetation52. Moreover, the changing climate and surface environment have significantly influenced the regional hydrologic cycle53, as evidenced by the spatiotemporal variations in runoff and ET54,55.

Cryosphere-hydrology model

The TP serves as a key area where the atmosphere, hydrosphere, biosphere, lithosphere, and cryosphere interact, with the hydrological cycle playing a central role in these interactions44,56,57,58. However, many existing studies exploring hydrological processes on the TP have been hindered by methodological constraints that do not fully consider these multi-sphere interactions. Considering these limitations, Wang et al.59 coupled a biosphere land surface model (SiB2) with a geomorphology-based hydrological model (GBHM) to develop the Water and Energy Budget-based Distributed Hydrological Model (WEB-DHM)59. This model considers water, energy, and CO2 fluxes within soil-vegetation-atmosphere systems at the grid scale, which not only incorporates energy conservation in land surface processes, but also delivers detailed representations of hydrophysical phenomena at high spatial resolutions59,60,61,62. The structure of the model at the grid scale is shown in Fig. S1. In addition, the WEB-DHM model has demonstrated its effectiveness in elucidating effects of cryosphere change on hydrological processes.

For snow simulation, the WEB-DHM incorporates an advanced enthalpy-based snow processes, employing the three-layer energy-balance snow parameterization from Simplified Simple Biosphere 3 (SSib3) and the prognostic albedo scheme from the Biosphere–Atmosphere Transfer Scheme (BATS)60,62,63,64,65. This module separates the snowpack into layers at each grid cell, with snow depths exceeding 5 cm being categorized into three layers; shallower depths, and those on canopy, are treated as a single layer (Fig. S2). This module considers the energy exchange between snow layers and the impact of solar radiation incidence angle and snow age on albedo. To sum up, this model could grant a meticulous depiction of snow physics, inclusive of phase transitions, compaction, albedo variations, temperature profiles, and meltwater runoff for each layer60,62,63,64,65.

Simultaneously, the model’s empirical frozen ground parameterization has been enhanced, encapsulating frozen ground dynamics through a hydrothermal transfer parameterization using the Johansen thermal conductivity framework66 (Fig. S2). By employing enthalpy rather than temperature in energy balance calculations, this approach reduces uncertainties associated with phase change latent heat release and improves model stability61,67,68. As a result, the model effectively simulates freeze-thaw cycles at fine scales and basin-scale hydrological processes, while appropriately accounting for subgrid topographic variations67.

The WEB-DHM model also integrates an energy balance-oriented glacier module to address both clean and debris-covered glaciers, thus augmenting the model’s capabilities for glacier processing69. For clean glaciers, a similar three-layer water- and energy-balance scheme is used, whereas a single-level scheme that accounts for debris influence on energy absorption is employed for debris-covered glaciers.

Water balance method

Given the challenge of directly measuring ET at the basin scale, WEB-DHM ET data were indirectly validated against the residual of precipitation (P) minus observed runoff (Robs) and terrestrial water storage changes measured by the Gravity Recovery and Climate Experiment (GRACE) satellite (ΔSGRACE) as the reference value42. However, this water balance approach dependent on the GRACE satellite is only applicable to very large regions, and cannot account for glacier-related changes in complex mountain regions (e.g., Ganges and Indus basins). Consequently, this method was only employed to evaluate the simulated ET for the Yellow, Yangtze, and Brahmaputra (upstream of the Nuxia Station) basins, where the glacier coverage is minimal and the basin area is large.

Performance metrics

A suite of metrics was utilized to assess the performance of the model and reliability of various products; these included the Nash coefficient70, Bias, Root Mean Square Error (RMSE), and Correlation Coefficient (CC)71. These metrics are computed as follows:

$${Nas}h=1-\frac{\mathop{\sum }\limits_{i=1}^{n}{\left({X}_{i}-{{ref}}_{i}\right)}^{2}}{\mathop{\sum }\limits_{i=1}^{n}{\left({{ref}}_{i}-\overline{{ref}}\right)}^{2}}$$
(1)
$${Bias}=\frac{\mathop{\sum }\limits_{i=1}^{n}\left({X}_{i}-{{ref}}_{i}\right)}{\mathop{\sum }\limits_{i=1}^{n}{{ref}}_{i}}\times 100 \% $$
(2)
$${RMSE}=\sqrt{\frac{\mathop{\sum }\limits_{i=1}^{n}{\left({X}_{i}-{{ref}}_{i}\right)}^{2}}{n}}$$
(3)
$${CC}=\frac{\mathop{\sum }\limits_{i=1}^{n}\left({X}_{i}-\bar{X}\right)\left({{ref}}_{i}-\overline{{ref}}\right)}{\sqrt{\mathop{\sum }\limits_{i=1}^{n}{\left({X}_{i}-\bar{X}\right)}^{2}}\sqrt{\mathop{\sum }\limits_{i=1}^{n}{\left({{ref}}_{i}-\overline{{ref}}\right)}^{2}}}$$
(4)

here, \({{\boldsymbol{ref}}}_{{\boldsymbol{i}}}\) represents the reference value at the i moment, such as observed flow or ET calculated using the water balance method; X is the target variable being evaluated, and n is the number of samples; \(\overline{{\boldsymbol{ref}}}\) and \(\bar{{\boldsymbol{X}}}\) represent the mean values of the reference sequence \({{\boldsymbol{ref}}}_{{\boldsymbol{i}}}\) and target sequence \({{\boldsymbol{X}}}_{{\boldsymbol{i}}}\), respectively. The closer the values of Nash and CC are to 1, the higher the reliability of the target variable. Conversely, the closer the values of RMSE and Bias are to 0, the higher the reliability of the target variable.

Datasets used in model input and evaluation

The WEB-DHM model incorporated in situ observations, reliable remote sensing, and reanalysis data as input. These datasets comprised observed daily discharge, gridded meteorological products, vegetation dynamic products, DEM, as well as detailed maps of land cover, soil, and glacier. Additionally, diverse satellite-based products such as land surface temperature and snow depth were utilized to validate the model’s multiple outputs. Specific datasets for each headwater model are listed in Table 2 and illustrated in Fig. 2.

Fig. 2
figure 2

Flowchart of this study.

Full size image

Data products for cross-comparison with TPRED

In recent decades, there has been a marked increase in runoff and ET products, each characterized by distinct input data, spatiotemporal resolutions, and modeling methodologies. While in situ measurements provide exceptional accuracy, remote sensing, model simulations, and data fusion techniques also contribute significantly, despite their inherent uncertainties23,25,29,72,73,74,75. These uncertainties stem from global discrepancies in gauge coverage, difficulties in accurately representing hydrophysical processes within remote sensing or machine learning frameworks, and the constraints about model performance or the quality of model inputs37,38,41,76,77,78. A detailed summary of limitations associated with numerous state-of-the-art runoff and ET products in the TP has been conducted; please refer to the Supplement “Part 1.2” for more information.

Therefore, it is crucial to meticulously evaluate the appropriateness of different data products for specific applications within the study region. This study presents a comprehensive comparison of the TPRED dataset with current gridded runoff and ET datasets, including runoff products from ERA5, ERA5-land, CNRD V1.0, CRUN, TerraClimate, and ET products from ERA5, ERA5-land, REA, ET Monitor, and GLEAM, as detailed in Table 3. Table 3 highlights the recurring deficiencies in these existing datasets, such as insufficient observations within the TP, relatively coarse spatial resolution, and frequent omission of critical cryosphere-hydrophysical processes.

Table 3 Datasets or products for cross-comparisons with TPRED.
Full size table

Data Records

The TPRED has been uploaded to the Zenodo platform and is publicly accessible via https://zenodo.org/records/10060590 or https://doi.org/10.5281/zenodo.1006059079. The time range is 1998–2017 at a monthly scale, the spatial resolution is 5 km × 5 km, and the unit is mm/month. For the convenience of users, this data is stored in TIF format, with each combination of different watersheds forming a separate file. Each file contains two variables: either runoff or ET.

Technical Validation

Calibration and validation of WEB-DHM models

We conducted model calibration and validation utilizing both ground-based and satellite-based observations, with a primary focus on observed daily discharge (Fig. 2). Model calibration involves adjusting saturated hydraulic conductivity (Ksurface), soil anisotropy ratio (anik), two van Genuchten parameters (α and n) and so on. Among these, Ksurface and anik are identified as the most sensitive parameters influencing streamflow simulations within the WEB-DHM. Ksurface primarily modulates peak flow, whereas anik impacts both peak and recession flows60,66. The α and n adjust the soil hydraulic function, influencing water transport within the soil and consequently affecting the flood hydrograph. Furthermore, a sensitivity study for these model parameters in the WEB-DHM has already been provided (https://hess.copernicus.org/preprints/6/C3455/2010/hessd-6-C3455-2010-supplement.pdf)60.

Building on the calibration result, the optimized parameters are subsequently employed during the validation phase. Table 4 presents the metrics (e.g., Nash and Bias) for the model performance during the calibration and validation periods, demonstrating that the WEB-DHM produces satisfactory results across different headwaters.

Table 4 The model assessment during the calibration and validation periods for each of the study headwaters.
Full size table

Evaluation of the WEB-DHM outcomes against the observed reference

The rigorous assessment of outcomes of the WEB-DHM model was conducted across different basins. Figure 3 shows a comparison between the model-simulated monthly streamflow and the observed flow at the outlet of each headwater region, utilizing Nash and Bias scores as evaluation metrics. The WEB-DHM estimations of monthly river flow demonstrated satisfactory agreement with the observational data across seven headwaters. This is evidenced by Nash values spanning the range 0.77–0.93, and Bias values oscillating between −14.2% and 5.9%. The simulated discharge adeptly captured the temporal pattern of river flow, closely mirroring the peaks and troughs of observed streamflow in terms of magnitude and timing.

Fig. 3
figure 3

Normalized observed and simulated discharge at the outlet of each headwater region across the TP. The Z-score method is used to standardize the data series96.

Full size image

For the TPRED ET simulations, the model demonstrated satisfactory accuracy in monthly ET simulation for the Yellow, Yangtze, and Brahmaputra (above the Nuxia station) rivers, as gauged against ET calculations derived from the water balance method (Fig. 4). The Yellow River displayed the most robust performance, signified by the highest Nash value (0.77), and significantly lower Bias (9.0%) and RMSE (15.3 mm). Skill scores for both the Yangtze and Brahmaputra rivers were also within reasonable bounds. It is notable, however, that the skill scores for TPRED ET were marginally lower than those for runoff in these catchments, possibly owing to inherent errors in the water balance approach and the coarse spatial resolution of the GRACE data.

Fig. 4
figure 4

Comparison of the TPRED ET with the monthly reference ET. Reference ET was obtained from the water balance approach (=P – Robs – ΔSGRACE), averaged over the Yellow, Yangtze, and Brahmaputra headwater regions. The right bar signifies the density scale, where a higher value corresponds to an increased concentration.

Full size image

In summary, the credibility of WEB-DHM model was assessed for each headwater region by comparing the model’s outcomes with observed river discharge and ET references. The resulting high skill scores underscored the efficacy of the WEB-DHM model in accurately representing complex hydrological dynamics within cryosphere-influenced basins.

Inter-comparison between the TPRED and other datasets

To enhance the validation of the reconstructed TPRED, it was compared with various widely recognized datasets. The assessment was done through the multi-dimensional comparison, including the monthly and annual variability, as well as spatial patterns and magnitudes across different topography and elevations.

Verification in runoff and ET products on the monthly scale

The monthly characteristics of six runoff datasets were examined in comparison to observational benchmarks, as illustrated in Fig. 5. The monthly runoff distribution in all datasets exhibited a consistent seasonal pattern, with values peaking during summer and reaching nadirs through winter and spring (Fig. S3). The TPRED outperformed other runoff products in capturing the temporal dynamics of high and low runoff periods, especially in recognizing “double peaks” in the Yellow headwater occurring in July and October. Additionally, the TPRED achieved high skill scores in the majority of headwater regions, supported by CC values exceeding 0.94. The TPRED consistently showed the lowest RMSE values as well as the relatively small Bias values when compared against other datasets (Fig. 6). In addition to TPRED, the ERA5 exhibited the second-highest agreement with observed references; CRUN significantly underestimated peak values in the Ganges and Indus headwaters while overestimating values in other basins; TerraClimate demonstrated relatively poor performance in terms of TP.

Fig. 5
figure 5

Comparison of ___domain-averaged monthly runoff obtained from six products with observed references averaged over the period from 1998 to 2017. The study divides the observed discharge at the outlet of each headwater by the basin area to estimate the runoff reference.

Full size image
Fig. 6
figure 6

Performance metrics of each runoff product compared to observed reference. The specific indexes include (a) CC, (b) RMSE, and (c) Bias, computed based on the data presented in Fig. 5.

Full size image

Figure 7 illustrates the monthly distribution of ET derived from TPRED, in comparison with five alternative ET datasets. TPRED consistently captures the expected monthly ET peaks, though it exhibits slight underestimations in some rivers compared to other datasets. Notably, GLEAM ET values are generally lower than those from other products across most basins. The CC between TPRED and other ET products is exceptionally high on a monthly scale, with all basin values exceeding 0.935. This strong correlation indicates that TPRED’s monthly ET performance aligns well with anticipated standards (Fig. 8). The seasonal ET pattern follows the runoff pattern, with the highest ET occurring in summer and progressively lower values in autumn, spring, and winter (Fig. S4).

Fig. 7
figure 7

Comparison of ___domain-averaged monthly ET obtained from six products averaged over the period from 1998 to 2017.

Full size image
Fig. 8
figure 8

The CC between each ET product and TPRED on the month scale during 1998~2017.

Full size image

In brief, monthly comparisons revealed that the TPERD runoff corresponds superiorly to observed reference than other products in terms of magnitude and seasonality, and the TPRED ET estimates closely matched expected seasonal fluctuations.

Multi-dimensional comparison of runoff and ET products on the annual scale

Figure 9 compares the spatial distributions of the long-term mean annual runoff as depicted by TPRED and other products. The runoff patterns portrayed by TPRED are consistent with those observed in other datasets, showing significant spatial variability. Water-rich basins, such as the Ganges, southeast Brahmaputra, and southwest Indus, contrast sharply with less affluent regions, including the Yangtze, Yellow, and Mekong basins. Variability in runoff across these basins can be largely attributed to discrepancies in rainfall distribution, as precipitation remains the primary source of runoff replenishment. This relationship is evidenced by a significant positive correlation between precipitation and runoff, excluding glacier-covered areas (Fig. S5). Furthermore, TPRED demonstrates superior capability in depicting the spatial heterogeneity of water resources across complex landscapes, especially in areas with substantial variations in elevation.

Fig. 9
figure 9

Spatial distribution of annual runoff among six products averaged over the period from 1998 to 2017 (unit: mm).

Full size image

The annual average runoff from the TPRED and other datasets was compared and visualized multidimensionally, revealing concordance across basins in terms of magnitude, except for TerraClimate, which significantly underestimated runoff (Fig. 10a). TPRED runoff values varied from 140 mm (Yangtze) to 1011 mm (Ganges). Interannual anomalies showed similar trends among all datasets, with minima in 2006, 2009, and 2015, and common peaks in 2000 and 2003 (Fig. 10b). Probabilistic and altitudinal distributions demonstrated that, despite magnitude discrepancies in specific elevation ranges, the overall pattern of runoff variation with altitude was consistent across most products (Fig. 10c,d).

Fig. 10
figure 10

Multi-dimensional comparison on six runoff products at the annual scale across different headwater regions. (a) Mean annual runoff for the period 1998–2017. (b) Temporal variation of annual runoff anomalies, defined as the differences between annual runoff values and the multi-year average for each product. (c) Probability distribution. (d) Altitudinal distribution.

Full size image

A comprehensive comparison was performed on ET datasets, as illustrated in Figs. 11 and 12. The southern Ganges and western Indus headwaters exhibited the highest ET values, whereas the Yangtze showed the lowest annual average ET. The Yellow and Mekong presented higher ET values than expected, given their respective precipitation and runoff levels. In the contexts of spatial distribution, probabilistic behavior, and altitudinal gradients, the ERA5 dataset displayed the greatest congruence with TPRED. The GLEAM consistently underestimated the magnitude of ET and failed to accurately capture interannual variability. The ERA5 dataset demonstrated inverse behavior concerning interannual fluctuations compared to most other datasets. The ET Monitor’s performance was closely aligned with the REA, although it showed a slight overestimation. Generally, ET decreased with elevation, with TPRED indicating slight fluctuations between 3000 and 4000 meters.

Fig. 11
figure 11

Spatial distribution of annual ET among six products averaged over the period from 1998 to 2017 (unit: mm).

Full size image
Fig. 12
figure 12

Multi-dimensional comparison on six ET products at the annual scale across different headwater regions. (a-d) Similar to Fig. 11, but for ET.

Full size image

Collectively, the annual values of TPERD generated by the rigorously validated WEB-DHM model showed commendable alignment with other established products across multiple dimensions, including watershed-scale magnitudes, temporal deviations, and elevational gradients. TPERD, along with other products, demonstrated similar spatial variability in runoff and ET distributions. The superior spatial resolution of TPERD enabled a more detailed delineation of spatial heterogeneity in runoff and ET across topographically complex terrains compared to coarser-grained datasets.