the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 23 Jun 2026
Hierarchical Graph Networks for Seasonal Forecasts of Terrestrial Water Storage Anomalies
Abstract. Fresh water availability is critical for ecosystems, agriculture, industry, and human communities. Anticipating drought conditions benefits from forecasting changes in terrestrial water storage (TWS), the total water stored on land across all compartments, including groundwater, rivers, glaciers, and soil moisture. While individual compartments, such as groundwater, are difficult to observe directly at large scales, TWS integrates their combined changes and can be measured globally through satellite gravimetry. Since 2002, the GRACE and GRACE Follow-On (GRACE-FO) missions have delivered monthly, global estimates of terrestrial water storage anomalies (TWSA), deviations from a long-term mean, making TWSA the most accessible large-scale indicator of hydrological change. Predicting TWSA is nonetheless challenging as it reflects processes operating at vastly different temporal and spatial scales. We present HiGNN-LSTM, a hierarchical graph neural network that represents the Earth across two spatial scales, coupled with an LSTM module to forecast global TWSA for up to six months ahead. As a proof of concept, we show that the hierarchical graph neural network can automatically generate meaningful input features for TWSA forecasting from ERA5 climate variables, without requiring manual predictor selection or lag-correlation analysis. Trained on a GRACE-like reconstruction of TWSA extending from 1979 to 2020, the model substantially reduces the one-month-lead RMSE relative to a seasonal climatology baseline (1.83 cm vs 3.70 cm) and consistently outperforms a ConvLSTM across the full six-month horizon. Skill over climatology shrinks at longer leads and is lost by six months, indicating that most of the gain concentrates at short leads. Evaluation against GRACE- and GRACE-FO-derived TWSA highlights the difficulty of transferring a model trained on reconstructed TWSA to satellite-derived observations.
- Preprint
(6059 KB) - Metadata XML
- BibTeX
- EndNote
Status: final response (author comments only)
- RC1: 'Comment on egusphere-2026-2312', Anonymous Referee #1, 04 Aug 2026
-
RC2: 'Comment on egusphere-2026-2312', Anonymous Referee #2, 11 Sep 2026
General Assessment
The manuscript presents HiGNN-LSTM, which is a two-level graph neural network coupled with LSTM, trained on the Li et al. (2021) TWSA reconstruction with ERA5 forcing, to forecast global TWSA at lead times of one to six months. The writing is clear and well-organized, and code and processed data are public which are all valuable.
My main concern is that the central claim of this study mentioned in abstract, introduction, and conclusion that the GNN automatically learns useful predictors from ERA5 is not tested by the experiments presented. In addition, several details of the model structure, implementation and workflow and hydrological interpretations of the results are largely absent. Therefore, I recommend major revision and below are my major concerns:
Major concerns
1. The model isn’t tested by isolating it from ERA5. Is there any evidence that adding ERA5 substantially improves the model? Retrain the model with only TWSA and month of the year and see if the performance drops. Then check if the ERA5 variables have contribution to the model or just adding complexity to the model. If the paper's claim is to improve performance by ERA5, there should be a clear comparison with and without ERA5.
2. The evaluation against GRACE observations (Table 3, Sect. 4.3) lacks sufficient baselines and diagnostics.
2a) Inconsistency between reconstruction and GRACE evaluation
Table 2 and Table 3 need meaningful and strong interpretation based on the comparison made against observation. The model loses to climatology at lead 6 on the reconstruction (3.77 vs. 3.69 cm) but beats climatology at every lead against GRACE, including lead 6 (6.96 vs. 7.96 cm). It is not clear why or how the model performs worse than the baseline on the reconstruction but better on the GRACE observation. It should be clearly stated how Table 3 climatology was computed.
From the released code in the repository evaluation_HiGNN.py, it seems that the climatology is computed from the reconstruction training period get_climatology() and reused unchanged for the GRACE evaluation. If so, then Table 3 compares the model against climatology derived from different smoother product which would explain the values in the table.
2b) The error in Table 3 stops growing as the lead time increases, and this is never discussed
Table 2 almost doubles across the horizon (1.83 → 3.77 cm), as expected. Table 3 goes 5.59 → 6.63 → 6.98 → 6.98 → 6.89 → 6.96 cm: it rises for two months, then flattens and even dips at lead 5. This suggests a fixed error floor near 7 cm caused by the mismatch between training data and observations rather than by forecast difficulty. Please discuss this, and clarify whether leads 3 to 6 in Table 3 measure forecast skill or a lead-independent error floor.
2c) The ConvLSTM baseline is missing from the GRACE evaluation
Table 2 compares HiGNN-LSTM with a ConvLSTM and climatology but in Table 3, ConvLSTM is dropped without any explanation. Since the ConvLSTM supports the claim at lines 214–215 that GNN features beat local convolutions, its absence removes that comparison from the one test of whether the advantage holds on real observations.
3. No uncertainty quantification is reported.
It appears from table 2 and table 3 that results are reported based on a single run without any seed variability, confidence interval, and significance test. Appendix B (l. 334–338) states 0.05 to 0.12 cm difference “confirming” while treating differences of at most 0.05 cm as negligible, which is self-contradictory with no variance estimate.
4. The overall workflow what is described in Figure 1 is ambiguous
It is not clear where the autoregressive forecast and the TWSA input sequence are implemented in the workflow, and what part of the model contributes to the autoregressive component. If it gets fixed for 12 months then it is not autoregressive. It should be formulated clearly if it is fixed offset prediction or recursive prediction from previous month. The figure itself with a minimal formulation of the method is ambiguous. There is only one equation in the manuscript which has no role in the model description. Although Figure 1 is useful for the overall workflow and spatial pipeline but does not have enough information about the temporal formulation. I would ask the authors to correct or explain the arrow (autoregressive label) and to add a short algorithm block, or a small set of equations, making the forecast procedure explicit. In particular, what the LSTM receives at each lead, whether the 12-month history is held fixed, what the residual is added to, and whether the recurrent state is carried across leads.
5. The role of seasonality in model performance
The seasonal cycle is kept in target, and the month of the year is given as an input to the model. This model will then reproduce seasonality very well and because of this the correlation will be very high in the basins which are dominated by the seasonal cycles like Amazon (r=0.99, Fig. 4). In these cases, seasonal cycle accounts for most of the variance. This is mentioned in manuscript section 4.2 but the model’s skill on anomaly components hasn’t been reported. Please report anomaly correlation with the seasonal cycle removed from both forecast and target, using a training-period climatology, per basin and per lead, together with skill scores against climatology and persistence.
Closing assessment
I recommend major revision. The main issues are the need to provide evidence for the model’s usefulness and performance in comparison to the baseline, uncertainty quantification of the model, clear description of the model and workflow in different mentioned parts. If these issues are addressed, the manuscript could make a valuable contribution to literature.
Citation: https://doi.org/10.5194/egusphere-2026-2312-RC2
Data sets
HiGNN_LSTM: Dataset for training and evaluating TWSA forecasts Viola Steidl https://doi.org/10.5281/zenodo.19664592
Model code and software
HiGNN-LSTM Viola Steidl https://github.com/viola1593/HiGNN-LSTM
Viewed
| HTML | XML | Total | BibTeX | EndNote | |
|---|---|---|---|---|---|
| 99 | 78 | 10 | 187 | 12 | 10 |
- HTML: 99
- PDF: 78
- XML: 10
- Total: 187
- BibTeX: 12
- EndNote: 10
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
I appreciate the substantial technical effort involved in developing and implementing this global hierarchical graph model. However, I have more fundamental concerns about the study's hydrological framing. The manuscript is primarily organized around a machine-learning architecture and then applies that architecture to TWSA, but it does not yet demonstrate that the resulting framework addresses a clearly formulated hydrological question. In particular, the central interpretation that the graph learns useful hydrological teleconnections is not supported by the current experiments. At present, the manuscript demonstrates that the proposed network can reproduce a reconstructed TWSA product more clearly than it demonstrates a new hydrological understanding or a useful advance in seasonal hydrological prediction. I therefore consider major revision necessary, and I am uncertain whether the required changes can be achieved without a substantial redesign of the study.
1) The current interpretability analysis does not yet establish which hydrological information the model has learned. The attribution analysis only shows that the output is sensitive to the GNN latent features, but it does not identify where the predictive information comes from. TWSA is primarily controlled by local water-balance processes and storage memory, and remote climate connections may provide additional information in specific regions and seasons. I suggest adding process-oriented attribution experiments and determining whether the learned remote connections correspond to known hydroclimatic relationships.
2) The practical hydrological value of forecasting total TWSA over a six-month horizon needs to be clarified. A technically accurate prediction of this integrated signal does not necessarily translate into useful forecasts of drought, groundwater deficit, reservoir conditions, or other hydrological impacts. I suggest that the authors define the intended hydrological application and evaluate metrics relevant to that application.
3) More importantly, because the seasonal cycle is retained in the target and the month of the year is explicitly provided as an input, the reported RMSE and correlation can be strongly influenced by the model’s ability to reproduce seasonality. For seasonal hydrological prediction, skill in the anomaly component is usually more informative than skill in the full signal. I recommend evaluating deseasonalized TWSA anomalies directly. The results should be reported separately for each lead time and basin.
4) The relationship between the predictor data and the reconstructed TWSA target requires a more careful discussion. The Li et al. (2021) reconstruction was itself produced using climatic and hydrological predictors and machine-learning methods, but the present model uses ERA5 climatic and hydrological variables to predict that reconstruction. Good performance on the reconstructed product may partly reflect HiGNN-LSTM's ability to reproduce relationships already built into that product. Ideally, the author may wish to repeat the analysis with more than one reconstruction product or provide a stronger observation-based evaluation.
5) Right now, it seems every ocean cluster is connected to every basin. I cannot immediately understand this strong structural assumption of universal ocean-to-basin connectivity. I suggest testing graph variants without ocean nodes, with only local basin connections, and with river-network connections. Otherwise, I don't think improved performance can be interpreted as evidence that the model has identified physically meaningful teleconnections.
Specific comments:
L24: However, these effects are not explicitly represented in the model inputs or graph structure, and their treatment in the reconstructed TWSA target remains unclear. The authors should clarify the extent to which the reconstruction preserves signals of human interventions. They should also discuss how detrending may remove persistent anthropogenic storage changes and how this affects transferability to GRACE observations and performance in strongly managed basins.
L60: Please specify which parts are hydrologically informed. The basin aggregation might be hydrologically motivated, but the nearest-neighbor and universal ocean-to-basin connections appear primarily geometric. Not sure of the value of basin aggregation here.
L69: For a hydrology journal, the manuscript should formulate one or more hydrological hypotheses in addition to this architectural objective. Otherwise, I have to say this manuscript is more suitable for a general Earth data science journal.
L82: Please explain why remote relationships are expected to contribute materially to TWSA prediction relative to local storage memory and local water-balance inputs with literature.
L92: Please specify the basin-selection procedure used to obtain the basin nodes. How are grid cells that intersect multiple basins treated?
L97: Connecting each basin to its geographically nearest neighbors does not necessarily represent hydrological interaction. Please explain how information transfer between hydrologically unrelated neighboring basins should be interpreted.
L133: Because the Li et al. (2021) reconstruction uses climatic and hydrological predictors, please compare its predictor set with the variables used in this study.
L145: Removing linear trends changes the scientific target. Please explain which hydrological signals are removed by detrending; some regions may be affected by long-term groundwater depletion or glacier mass loss.
L161: The CSR mascon product is just distributed on a 0.25-degree grid, but this grid spacing does not represent its independent spatial resolution. Please make sure this sentence is precise.
Section 4.1: Aggregating attribution values over all nodes removes the spatial and temporal information that is most relevant to the graph hypothesis. Please provide basin-specific examples and show how attribution changes across input months and lead times.
L191: All four latent dimensions having nonzero attribution by themselves does not demonstrate that they contain distinct or useful information. A latent-dimension removal experiment would provide more direct evidence.
L261: I cannot understand here. Predicting the full signal may be harder in absolute RMSE terms, but the strong seasonal cycle can make correlation and normalized skill easier to obtain. The comparison should use common anomaly-based metrics.
L290: The poor performance in anomaly-dominated basins may also reflect human water management, groundwater abstraction, snow processes, or reconstruction errors. The discussion currently focuses mainly on seasonality and graph connectivity.
L306: The experiments show that the full model outperforms ConvLSTM and an embedder-only variant, but they do not identify the learned connections as teleconnections.