the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 10 Aug 2026
A Gauge-Based Monthly Natural Discharge Reconstruction Dataset for the Eurasian Arctic, 1952–2025
Abstract. River discharge provides one of the clearest measures of how Arctic freshwater systems respond to climate change, yet long-term and continuous monitoring remains challenging. Here we present a gauge-based monthly natural discharge reconstruction dataset for 94 gauging stations in the Eurasian Arctic from 1952 to 2025. The dataset was generated using a basin-scale framework driven by precipitation and air temperature and aided by snow information. Terrestrial water storage (TWS) serves as the key link in the framework. Historical TWS was first reconstructed using an improved state-update empirical model constrained by TWS observations from the Gravity Recovery and Climate Experiment (GRACE), and monthly natural discharge was then estimated from reconstructed TWS using a refined runoff–TWS relationship calibrated against gauge discharge observations. The reconstructed monthly natural discharge agreed well with observations at minimally regulated gauges and during pre-regulation periods at regulated gauges, with median monthly Kling–Gupta efficiency (KGE) values exceeding 0.88. When monthly values were aggregated to annual discharge, the reconstruction still captured observed variability at minimally regulated gauges, with a median annual KGE of 0.73. The dataset generally agreed better with gauge observations than discharge estimates derived from existing GRACE-like TWS products and an existing gridded runoff reconstruction product. Example applications illustrate spatially heterogeneous trends in annual discharge and seasonal allocation over the past seven decades and show how the dataset can be used to assess the effects of human regulation on river discharge. This dataset provides a long-term natural-discharge baseline for characterizing historical discharge variability and assessing the roles of climate variability and human regulation in poorly monitored Arctic basins. The reconstructed monthly natural discharge dataset is available at https://doi.org/10.5281/zenodo.21157816 (Liu et al., 2026).
- Preprint
(5828 KB) - Metadata XML
-
Supplement
(2229 KB) - BibTeX
- EndNote
Status: open (until 25 Sep 2026)
-
RC1: 'Comment on essd-2026-531', Anonymous Referee #1, 10 Sep 2026
reply
-
AC1: 'Reply on RC1', Bingshi Liu, 13 Sep 2026
reply
==================================================
This comment is a preliminary response submitted during the public discussion to clarify several points raised by the reviewer and to outline our planned revisions. If we are invited to submit a revised manuscript after the discussion, we will revise the manuscript and the deposited dataset files in light of all reviewer comments and provide additional detailed responses to points that require further clarification or are addressed through manuscript and dataset revisions.
===================================================We thank the reviewer for the detailed and constructive assessment of our manuscript and deposited dataset.
First, we would like to clarify that the reconstructed natural discharge has already been evaluated against gauge discharge observations through multiple validation experiments in the current manuscript. These include independent testing at minimally regulated gauges, where the evaluation periods are outside the parameter-identification periods, and pre-regulation validation at regulated gauges, both presented in Section 4.1. In Section 4.2, we further evaluated the reconstruction at the annual-discharge scale and using apportionment entropy, which provides information on year-to-year variations in discharge magnitude and seasonal allocation. Therefore, the current evaluation is not limited to the mean seasonal cycle.We nevertheless agree with the reviewer that the validation design of the reconstruction framework can be made clearer. In the revised manuscript, we will add separate evaluations of the seasonal and non-seasonal signals in reconstructed monthly natural discharge against gauge discharge observations for both the training periods and the independent testing periods. The non-seasonal signal will be obtained by removing the mean seasonal cycle. These analyses will further assess whether the reconstruction can reproduce variability beyond the mean seasonal cycle, including interannual anomalies, long-term variability, and low-flow conditions.
Q1: Scientific significance and added value (Introduction, Sections 4–5, and Conclusions). The geographical coverage and length of the reconstruction are potentially useful, but they do not by themselves establish its scientific contribution. How much information does the intermediate TWS reconstruction add? Please compare the framework with a monthly climatology and a parsimonious discharge model using the same meteorological inputs and calibration periods. Evaluate these comparisons on independent observations, especially for interannual anomalies, low flows, and trends. This would establish whether the additional empirical structure provides meaningful predictive information and clarify which applications the dataset can support.
R1: We agree that geographical coverage and record length alone are insufficient to establish the scientific contribution of the dataset. The purpose of this dataset is to provide a model-based natural-discharge baseline that complements observed discharge records. Observed discharge represents actual river flow under the combined influence of climate variability and direct human regulation. In contrast, the reconstructed dataset estimates the discharge expected from climatic and snow-storage conditions without explicit human-regulation inputs. It can therefore help characterize discharge changes driven mainly by climatic and snow-related factors. Such a baseline is important for studying long-term freshwater variability and for examining links, interactions, and feedbacks among climate, hydrological processes, and Arctic ecosystems.
We would also like to clarify that the current manuscript does not evaluate only the agreement in the mean seasonal cycle. Section 4.1 evaluates reconstructed monthly natural discharge against gauge discharge observations at minimally regulated gauges during independent testing periods outside parameter identification and during pre-regulation periods at regulated gauges. Section 4.2 further evaluates annual discharge and apportionment entropy, which provide diagnostics of year-to-year variability in discharge magnitude and seasonal allocation. Therefore, the current evaluation already includes information beyond the mean seasonal cycle.
We understand the reviewer’s reference to monthly climatology as a concern that the reconstruction skill may be dominated by the strong seasonal cycle of Arctic rivers. To address this point more clearly, we will revise the validation section to explicitly separate the evaluation of seasonal and non-seasonal signals. The non-seasonal signal will be obtained by removing the mean seasonal cycle, and the reconstructed and observed non-seasonal components will be compared during both the training periods and the independent testing periods. This additional analysis will more directly show whether the reconstruction captures variability beyond the mean seasonal cycle.
Regarding the reviewer’s suggestion to compare the framework with a parsimonious discharge model, we understand that the purpose is to examine whether the proposed framework provides information beyond a simple meteorological reconstruction. We respectfully note that developing, calibrating, and optimizing an additional precipitation-temperature discharge model for all 94 gauges would introduce an additional modelling exercise that is not the primary objective of this dataset paper. The main objective of this study is to develop, document, and evaluate a TWS-based framework for reconstructing historical natural discharge.
Nevertheless, we agree that comparison with meteorology-based reconstructions is useful. We would like to clarify that the current manuscript already includes an external benchmark against GRUN (Ghiggi et al., 2019) in Section 5.2. GRUN is an independent gridded runoff reconstruction product driven mainly by precipitation and air temperature, and it does not explicitly represent direct human regulation or glacier-mass changes in the way considered here. It therefore provides a useful external reference for comparing our gauge-based reconstruction with an existing meteorology-based runoff reconstruction product.
Ghiggi, G., Humphrey, V., Seneviratne, S. I., & Gudmundsson, L. (2019). GRUN: an observation-based global gridded runoff dataset from 1902 to 2014. Earth System Science Data, 11(4), 1655–1674. https://doi.org/10.5194/essd-11-1655-2019Q2: Temporal coverage and changing input availability (Tables 1–3, Section 3.3, and deposited files). Reconstructed discharge and uncertainty are missing for January 1952 at all 94 gauges in both files; complete reconstructed coverage therefore begins in February 1952. Please supply these estimates or correct the stated coverage and explain the initialization procedure. Several input products end in 2014, 2017, 2019, 2023, or 2024, CPC begins in 1979, and neither listed snow-cover product extends through 2025. Explain how the stated fixed ensemble was generated throughout 1952–2025.
R2: We thank the reviewer for identifying the missing January 1952 values. In the runoff–TWS relationship, discharge estimation requires terrestrial water storage (TWS) information for both the current month and the preceding month. Because the effective reconstructed TWS series begins in January 1952, the January 1952 discharge estimate lacks the required December 1951 TWS state. Should a revision be invited, we will either reprocess the dataset by initializing the reconstructed TWS series from an earlier month so that the January 1952 discharge estimates can be derived consistently, or revise the manuscript, README, and file attributes to state clearly that the reconstructed discharge records start in February 1952.
We agree that the availability of input products and the ensemble construction should be described more clearly. Each forcing and auxiliary-data combination can be reconstructed only over the period common to its component datasets. Therefore, the final ensemble-mean product was generated by averaging the reconstruction members available in each month. This allows the ensemble-mean discharge product to cover 1952–2025, although the number and identity of ensemble members may vary through time. Because all members are constrained by the same observational datasets, they are expected to be broadly comparable. Nevertheless, we will revise the manuscript to avoid implying that a fixed set of ensemble members is available for the entire period, and we will examine and present member consistency during overlapping periods in the supplementary.
In addition, the time coverage of the ERA5-Land snow-cover product in Table 3 was incorrectly stated because of a typographical error. The correct coverage is 1950–2025.Q3: Physical meaning of the storage components (Section 3.1 and Eqs. (1)–(4)). Please define physically what is included in “solid TWS” and “liquid TWS”. Solid water storage can include snow, glacier ice, and ground ice, whereas the implementation approximates it by SWE. Consequently, subtracting SWE from GRACE TWS does not necessarily isolate physically liquid water storage, particularly in permafrost regions. Explain the treatment of frozen soil water, ground ice, glaciers, and surface-water ice, and distinguish operational model variables from actual hydrological storage compartments. The storage-offset parameter dw in Eq. (1) is discussed later, but should be defined immediately. Does it represent an effective fitted baseline, and what prevents TWSliquid+dw from becoming physically inadmissible? A fitted offset should not be interpreted as an independently established absolute storage volume. Eq. (2) estimates net monthly SWE depletion, which is not equivalent to snowmelt when snowfall and sublimation occur simultaneously. Explain this approximation and its consistency with the additional solid-precipitation contribution in Eq. (1).
R3: We thank the reviewer for this comment. We agree that the descriptions of “solid TWS” and “liquid TWS” in the current manuscript should be clarified further before introducing the equations, and that the rationale for approximating solid TWS by snow water equivalent (SWE) should be explained more explicitly.
Conceptually, solid TWS may include snow, glacier ice, ground ice, and other frozen water components. In this study, we approximate solid TWS with SWE because long-term, spatially continuous SWE datasets are available, whereas comparable datasets for other solid-water components remain limited. Previous studies have also shown that seasonal snow represents the dominant active solid-water component in non-glaciated northern high-latitude basins and contributes strongly to seasonal TWS variability (Trautmann et al., 2018). Accordingly, in this framework, liquid TWS denotes the remaining TWS component after removing the SWE-based solid-storage proxy. We will clarify this rationale in Section 3.1. The limitation associated with using SWE to approximate solid TWS has already been discussed in Section 5.4 of the current manuscript.
We will define the storage-offset parameter (dw) immediately after Eq. (1). We will clarify that (dw) is an empirical fitted offset used to keep (TWSliquid + dw) positive in the model, and that it should not be interpreted as an independently determined absolute storage volume. We will also revise the description of Eq. (2) by referring to this term as positive net monthly SWE depletion rather than gross snowmelt. The solid-precipitation term in Eq. (1) will be clarified as an empirical contribution of current-month solid precipitation to liquid-storage changes, which may also partly absorb uncertainties related to precipitation-phase partitioning.
Trautmann, T., Koirala, S., Carvalhais, N., Eicker, A., Fink, M., Niemann, C., & Jung, M. (2018). Understanding terrestrial water storage variations in northern latitudes across scales. Hydrology and Earth System Sciences, 22(7), 4061–4082. https://doi.org/10.5194/hess-22-4061-2018.Q4: Interpretation as natural discharge (Sections 3.2–3.3, 4, and 5.4). Calibration against pre-regulation discharge does not ensure that the entire framework is free from regulation effects. The TWS model is calibrated against GRACE observations from 2003–2012, which may contain reservoir-storage signals. Explain how these signals are excluded or how their influence on the fitted parameters is assessed. Similarly, the statement in lines 169–171 that non-target signals are absent from reconstructed TWS requires justification. The distinction between a climate-driven statistical reconstruction and a natural-discharge counterfactual should be made explicit. Fixed empirical relationships may not preserve changes caused by evolving permafrost conditions, subsurface connectivity, or lake storage. Reconcile the presentation of Salekhard and Kyusyur as minimally regulated examples in Fig. 4c with their treatment as regulated outlets in Section 6.2.
R4: We thank the reviewer for this comment. We will revise the definition of “reconstructed natural discharge” in line 85 of the current manuscript. In the revised manuscript, reconstructed natural discharge will be described as a model-based, climate-driven natural-discharge baseline, that is, the discharge expected in the absence of direct human regulation. It should not be interpreted as a fully process-resolved counterfactual simulation that removes all effects of human activities and landscape changes.
We agree that GRACE TWS observations may contain reservoir-storage signals. In the Eurasian Arctic, such signals during 2003–2012 are mainly associated with several large reservoirs in the upper Yenisei basin and the Vilyuy tributary of the Lena basin, and may influence the GRACE TWS observations for contributing areas corresponding to up to about 10 of the 94 gauges in this study. The original analysis did not explicitly remove these signals. In the revision, we will use available data sources to assess the relative contribution of large-reservoir storage variations to GRACE TWS changes in the corresponding contributing areas and evaluate their potential influence on TWS parameter identification. We will also revise the statement that reconstructed TWS contains no non-target signals. The revised manuscript will state that direct human regulation is not explicitly represented as an independent driver because we generate reconstructed TWS from precipitation, air temperature, and snow-related inputs.
Regarding permafrost-related hydrological changes, this limitation has already been discussed in Section 5.4 of the current manuscript, including possible effects of thermokarst lake changes in northern lowlands. We will make this limitation clearer.
We also thank the reviewer for pointing out the ambiguity related to Salekhard and Kyusyur. The “minimally regulated” classification used in Figure 4c refers to gauge-level evaluation, where observed discharge shows limited or weakly detectable regulation effects and is therefore suitable for assessing reconstructed natural discharge. This classification does not imply the complete absence of reservoirs or other human interventions within the entire upstream basin. Section 6.2 was intended as an example application at major river outlets, not as a reclassification of these gauges as regulated. We will revise the section title and wording, for example to “Potential regulation-related departures at major river outlets”, and avoid treating all outlet gauges analyzed in this section as regulated gauges.Q5: Monthly versus annual performance and validation independence (Section 4 and Figs. 3–6). Please explain why median KGE decreases from approximately 0.88 at the monthly testing scale to 0.73 at the annual scale. The comparison should also use consistent independent periods: monthly testing excludes calibration observations, whereas the annual and entropy evaluations include them. My calculation reproduced the annual statistics at selected gauges, but Khatanga’s near-perfect correlation in Fig. 5 is based on only three complete hydrological years, 1965–1967. This result should be explicitly flagged and should not be presented as strong evidence of interannual reconstruction skill.
R5: We thank the reviewer for this comment. The decrease in median KGE from the monthly testing scale to the annual scale is expected because monthly discharge is strongly controlled by the seasonal cycle. After monthly values are aggregated to annual discharge, the dominant seasonal cycle is removed, and the evaluation depends more strongly on interannual variability. In addition, the number of annual samples is much smaller than the number of monthly samples, making annual KGE more sensitive to individual anomalous years and to the correlation, variability, and bias components of the metric. We will add this explanation to the revised manuscript.
We also agree that annual discharge and apportionment entropy statistics that include calibration years should not be described as fully independent validation. In the current manuscript, all available complete hydrological years were used to maximize the use of sparse gauge observations. In the revision, we will separate independent testing-period evaluations from full-period descriptive evaluations where data availability permits.
We thank the reviewer for noting the limited number of complete hydrological years at Khatanga. We agree that a near-perfect correlation based on only three complete hydrological years should not be interpreted as strong evidence of interannual reconstruction skill. We will remove or replace Khatanga as an illustrative example in the revised figure.Q6: Uncertainty estimation (Section 3.3 and Eq. (6)). Define N and the ensemble indexing. The N−2 term corresponds to the variance of an ensemble mean under an independence assumption, which requires justification because members share observations, forcing information, and model structure. Clarify whether the intended quantity is uncertainty in the ensemble mean or predictive uncertainty in reconstructed discharge. Averaging the 24 TWS reconstructions before discharge estimation does not explicitly propagate their spread through the nonlinear discharge model. Propagate uncertainty through both modelling stages, address dependencies and structural error, and evaluate uncertainty-interval coverage against withheld observations.
R6: We thank the reviewer for identifying the need for clearer uncertainty definitions. In the current manuscript, Eq. (6) was intended to describe the uncertainty of the ensemble mean of reconstructed natural discharge within the adopted ensemble framework. We will revise the text to define N and the ensemble indexing explicitly. Here, N denotes the number of ensemble members used in the averaging, which is 10 for the reconstructed discharge product in this study.
We agree that the independence assumption among ensemble members should not be overinterpreted, because members share observational constraints, forcing information, and model structure. We will therefore clarify that the reported uncertainty mainly reflects the spread of the adopted reconstruction ensemble and the uncertainty of its ensemble mean, rather than a complete predictive uncertainty including all possible structural errors. We will also clarify that uncertainty from the reconstructed TWS is partly propagated into discharge reconstruction, because each SWE product corresponds to a TWS reconstruction based on the same SWE product, and differences among these TWS reconstructions are carried forward into the discharge ensemble.Q7: Methodological reproducibility and TWS evaluation (Sections 2–3). Essential details are missing: the MCMC likelihood, priors and parameter bounds, convergence assessment, normalization, initial conditions and spin-up, spatial resampling, and basin aggregation. Because reconstructed TWS is the central link in the framework, evaluate it directly against GRACE observations outside 2003–2012, including seasonal and interannual components. Good discharge agreement alone does not establish that this intermediate state is represented reliably; fitted discharge coefficients could compensate for errors in reconstructed storage.
R7: We thank the reviewer for this constructive comment. In the revised manuscript and Supplement, we will add methodological details.
In the current manuscript, we did not present the TWS reconstruction results in detail because the final dataset and main evaluation focus on reconstructed discharge. However, we agree that TWS is the key link in the framework and that we should evaluate its reliability directly. Therefore, in the revised manuscript, we will add an evaluation of reconstructed TWS against GRACE and GRACE-FO TWS observations outside the 2003–2012 parameter-identification period.
We agree that fitted discharge coefficients may partly compensate for errors in reconstructed storage. However, the comparison in Section 5.1 also shows that not all discharge estimates derived from existing GRACE-like TWS products agree well with gauge discharge observations. This indicates that fitted discharge coefficients cannot fully mask errors or differences in TWS reconstructions.Q8: Dataset documentation and usability (Section 7 and deposited files). The files lack machine-readable calibration periods, regulation flags, validation statistics, and ensemble membership information. Add these fields and a README describing processing and quality limitations. The NetCDF time variable is a two-row year/month array with units “year month”, rather than a conventional calendar-aware time coordinate. Catchment boundaries are available only in MAT format, and MATLAB string objects are not directly decoded by SciPy’s standard reader.
R8: We thank the reviewer for these useful suggestions on dataset documentation and usability. We agree that the deposited files should be more self-documented and easier to reuse. In the revision, we will improve the README and file metadata by adding key information needed to interpret the dataset.
We will also revise the NetCDF file to use a more conventional calendar-aware time coordinate. Where feasible, we will provide or document additional supporting information needed to reproduce the key analyses in the manuscript.Q9: Trend inference and attribution (Section 6). Explain how serial dependence, multiple testing across gauges, and reconstruction uncertainty are considered in trend significance. Demonstrate whether observed trends at minimally regulated gauges are reproduced over common periods. The mapped reconstruction trends may primarily reflect the forcing products and assumed stationary relationships. Observed-minus-reconstructed differences contain model error as well as possible regulation effects. In Figs. 11–12, regulated-gauge anomalies are standardized using training-period variability, whereas background anomalies use full-period variability. Use comparable independent reference periods and assess sensitivity before attributing departures to regulation.
R9: We thank the reviewer for this comment. We agree that trend significance should consider serial dependence, multiple testing among gauges, and reconstruction uncertainty. In the revised manuscript, trends will be estimated using Sen’s slope, and trend uncertainty will be assessed using a moving-block bootstrap to partly account for serial dependence. We will also clarify that the analyses in Section 6 are intended as example applications of the reconstructed natural-discharge dataset rather than as a formal attribution study. Gauge-level significance will therefore be interpreted conservatively, and observed-minus-reconstructed differences will be described as possible regulation-related departures rather than direct estimates of anthropogenic regulation effects.
To provide an observational diagnostic of trend reliability, we will compare observed and reconstructed annual-discharge trends over common observation periods at minimally regulated gauges with sufficient records. This analysis will assess whether the reconstruction reproduces observed trend behavior where gauge observations are available.
Regarding Figures 11–12, we thank the reviewer for pointing out the inconsistency in the standardization procedure. In the original analysis, anomalies at the four outlet gauges were standardized using standard deviations from their training periods because the remaining periods were affected by reservoir regulation. Using full-period standard deviations at these regulated gauges would include regulation-related departures in the reference variability and could inflate the denominator, thereby masking post-regulation signals. However, we agree that using full-period standard deviations for minimally regulated gauges and training-period standard deviations for regulated outlet gauges made the background range and outlet anomalies not directly comparable.
In the revised manuscript, we will recalculate Figures 11–12 using a consistent standardization procedure. Anomalies at both minimally regulated gauges and the four outlet gauges will be standardized by the standard deviation of observed-minus-reconstructed anomalies during the gauge-specific training period. This will ensure that the background range and outlet anomalies are expressed relative to the same type of baseline variability.Q10: Glacier representation (Fig. 1 and Section 3.1). Glacier extent appears visually exaggerated in Fig. 1. What is the source, reference date, and spatial resolution of the glacier boundaries? Clarify whether the white marks represent actual glacier polygons or enlarged location symbols. If they are symbols, state explicitly that they do not depict glacier area. Because limited glacier coverage is used to justify approximating solid TWS by snow storage, provide quantitative glacier fractions for the relevant catchments.
R10: We thank the reviewer for pointing out the potential ambiguity in Figure 1. The white marks in Figure 1 are enlarged glacier-location symbols derived from the Randolph Glacier Inventory (RGI) version 6.0 (Pfeffer et al., 2014) and do not depict actual glacier area. We will revise the figure caption to state this explicitly.
To support the statement of limited glacier coverage, we will report the fraction of 0.5° contributing-area grid cells intersecting RGI glacier polygons for each relevant catchment. This grid-cell-based metric is not intended to represent the actual glacier-covered area fraction. Because a grid cell is counted once it intersects any glacier polygon, it provides a conservative indicator of possible glacier influence at the 0.5° scale used in this study. If this fraction remains small, the actual glacier-covered area fraction is expected to be even smaller.Pfeffer, W.T., Arendt, A.A., Bliss, A., Bolch, T., Cogley, J.G., Gardner, A.S., Hagen, J.-O., Hock, R., Kaser, G., Kienholz, C., Miles, E.S., Moholdt, G., Mölg, N., Paul, F., Radić, V., Rastner, P., Raup, B.H., Rich, J., Sharp, M.J., The Randolph Consortium, 2014. The Randolph Glacier Inventory: a globally complete inventory of glaciers. J. Glaciol. 60, 537–552. https://doi.org/10.3189/2014JoG13J176.
We will also correct the technical errors listed by the reviewer, including typographical errors, inconsistent dataset-version descriptions, and figure-label mistakes.
Citation: https://doi.org/10.5194/essd-2026-531-AC1
-
AC1: 'Reply on RC1', Bingshi Liu, 13 Sep 2026
reply
Data sets
Monthly Natural Discharge Reconstruction for 94 Eurasian Arctic Gauges, 1952–2025 Bingshi Liu et al. https://doi.org/10.5281/zenodo.21157816
Viewed
| HTML | XML | Total | Supplement | BibTeX | EndNote | |
|---|---|---|---|---|---|---|
| 166 | 41 | 49 | 256 | 41 | 31 | 59 |
- HTML: 166
- PDF: 41
- XML: 49
- Total: 256
- Supplement: 41
- BibTeX: 31
- EndNote: 59
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
General comments
The manuscript presents a reconstruction of monthly natural discharge for 94 Eurasian Arctic gauges over 1952–2025. Long, consistent discharge records would be valuable for the region, particularly for investigating changes in annual runoff, seasonal distribution, and reservoir impacts. The manuscript is generally well organized, and providing both observations and reconstructed values facilitates evaluation.
I downloaded and inspected both files deposited on Zenodo. The numerical discharge and uncertainty arrays agree between formats after accounting for storage precision, and selected annual evaluation statistics are reproducible. However, these checks establish internal consistency rather than the reliability of the reconstructed historical changes. Important concerns remain regarding temporal coverage, uncertainty propagation, reproducibility, and the interpretation of the product as natural discharge.
I am also not yet convinced that the dataset’s scientific significance has been adequately demonstrated. Extending meteorological inputs through fitted storage–discharge relationships does not automatically provide reliable new information about historical river discharge. The authors need to demonstrate that the reconstruction captures independent hydrological variability beyond the seasonal cycle, adds value relative to other existing benchmarks, and remains informative outside its calibration conditions. These requirements are particularly important because the proposed applications concern long-term trends and attribution, whereas the framework does not explicitly represent several processes capable of changing those trends.
Specific comments
Scientific significance and added value (Introduction, Sections 4–5, and Conclusions). The geographical coverage and length of the reconstruction are potentially useful, but they do not by themselves establish its scientific contribution. How much information does the intermediate TWS reconstruction add? Please compare the framework with a monthly climatology and a parsimonious discharge model using the same meteorological inputs and calibration periods. Evaluate these comparisons on independent observations, especially for interannual anomalies, low flows, and trends. This would establish whether the additional empirical structure provides meaningful predictive information and clarify which applications the dataset can support.
Temporal coverage and changing input availability (Tables 1–3, Section 3.3, and deposited files). Reconstructed discharge and uncertainty are missing for January 1952 at all 94 gauges in both files; complete reconstructed coverage therefore begins in February 1952. Please supply these estimates or correct the stated coverage and explain the initialization procedure. Several input products end in 2014, 2017, 2019, 2023, or 2024, CPC begins in 1979, and neither listed snow-cover product extends through 2025. Explain how the stated fixed ensemble was generated throughout 1952–2025.
Physical meaning of the storage components (Section 3.1 and Eqs. (1)–(4)). Please define physically what is included in “solid TWS” and “liquid TWS”. Solid water storage can include snow, glacier ice, and ground ice, whereas the implementation approximates it by SWE. Consequently, subtracting SWE from GRACE TWS does not necessarily isolate physically liquid water storage, particularly in permafrost regions. Explain the treatment of frozen soil water, ground ice, glaciers, and surface-water ice, and distinguish operational model variables from actual hydrological storage compartments. The storage-offset parameter dw in Eq. (1) is discussed later, but should be defined immediately. Does it represent an effective fitted baseline, and what prevents TWSliquid+dw from becoming physically inadmissible? A fitted offset should not be interpreted as an independently established absolute storage volume. Eq. (2) estimates net monthly SWE depletion, which is not equivalent to snowmelt when snowfall and sublimation occur simultaneously. Explain this approximation and its consistency with the additional solid-precipitation contribution in Eq. (1).
Interpretation as natural discharge (Sections 3.2–3.3, 4, and 5.4). Calibration against pre-regulation discharge does not ensure that the entire framework is free from regulation effects. The TWS model is calibrated against GRACE observations from 2003–2012, which may contain reservoir-storage signals. Explain how these signals are excluded or how their influence on the fitted parameters is assessed. Similarly, the statement in lines 169–171 that non-target signals are absent from reconstructed TWS requires justification. The distinction between a climate-driven statistical reconstruction and a natural-discharge counterfactual should be made explicit. Fixed empirical relationships may not preserve changes caused by evolving permafrost conditions, subsurface connectivity, or lake storage. Reconcile the presentation of Salekhard and Kyusyur as minimally regulated examples in Fig. 4c with their treatment as regulated outlets in Section 6.2.
Monthly versus annual performance and validation independence (Section 4 and Figs. 3–6). Please explain why median KGE decreases from approximately 0.88 at the monthly testing scale to 0.73 at the annual scale. The comparison should also use consistent independent periods: monthly testing excludes calibration observations, whereas the annual and entropy evaluations include them. My calculation reproduced the annual statistics at selected gauges, but Khatanga’s near-perfect correlation in Fig. 5 is based on only three complete hydrological years, 1965–1967. This result should be explicitly flagged and should not be presented as strong evidence of interannual reconstruction skill.
Uncertainty estimation (Section 3.3 and Eq. (6)). Define N and the ensemble indexing. The N−2 term corresponds to the variance of an ensemble mean under an independence assumption, which requires justification because members share observations, forcing information, and model structure. Clarify whether the intended quantity is uncertainty in the ensemble mean or predictive uncertainty in reconstructed discharge. Averaging the 24 TWS reconstructions before discharge estimation does not explicitly propagate their spread through the nonlinear discharge model. Propagate uncertainty through both modelling stages, address dependencies and structural error, and evaluate uncertainty-interval coverage against withheld observations.
Methodological reproducibility and TWS evaluation (Sections 2–3). Essential details are missing: the MCMC likelihood, priors and parameter bounds, convergence assessment, normalization, initial conditions and spin-up, spatial resampling, and basin aggregation. Because reconstructed TWS is the central link in the framework, evaluate it directly against GRACE observations outside 2003–2012, including seasonal and interannual components. Good discharge agreement alone does not establish that this intermediate state is represented reliably; fitted discharge coefficients could compensate for errors in reconstructed storage.
Dataset documentation and usability (Section 7 and deposited files). The files lack machine-readable calibration periods, regulation flags, validation statistics, and ensemble membership information. Add these fields and a README describing processing and quality limitations. The NetCDF time variable is a two-row year/month array with units “year month”, rather than a conventional calendar-aware time coordinate. Catchment boundaries are available only in MAT format, and MATLAB string objects are not directly decoded by SciPy’s standard reader.
Trend inference and attribution (Section 6). Explain how serial dependence, multiple testing across gauges, and reconstruction uncertainty are considered in trend significance. Demonstrate whether observed trends at minimally regulated gauges are reproduced over common periods. The mapped reconstruction trends may primarily reflect the forcing products and assumed stationary relationships. Observed-minus-reconstructed differences contain model error as well as possible regulation effects. In Figs. 11–12, regulated-gauge anomalies are standardized using training-period variability, whereas background anomalies use full-period variability. Use comparable independent reference periods and assess sensitivity before attributing departures to regulation.
Glacier representation (Fig. 1 and Section 3.1). Glacier extent appears visually exaggerated in Fig. 1. What is the source, reference date, and spatial resolution of the glacier boundaries? Clarify whether the white marks represent actual glacier polygons or enlarged location symbols. If they are symbols, state explicitly that they do not depict glacier area. Because limited glacier coverage is used to justify approximating solid TWS by snow storage, provide quantitative glacier fractions for the relevant catchments.
Technical corrections