Machine-learning-based reconstruction of long-term global terrestrial water storage anomalies from observed, satellite and land-surface model data

Mandal, Nehar; Das, Prabal; Chanda, Kironmala

doi:https://doi.org/10.5194/essd-17-2575-2025

Articles | Volume 17, issue 6

https://doi.org/10.5194/essd-17-2575-2025

© Author(s) 2025. This work is distributed under
the Creative Commons Attribution 4.0 License.

https://doi.org/10.5194/essd-17-2575-2025

© Author(s) 2025. This work is distributed under
the Creative Commons Attribution 4.0 License.

Articles | Volume 17, issue 6

Data description paper

|

13 Jun 2025

Data description paper |

| 13 Jun 2025

Machine-learning-based reconstruction of long-term global terrestrial water storage anomalies from observed, satellite and land-surface model data

Nehar Mandal, Prabal Das, and Kironmala Chanda

Download

Final revised paper (published on 13 Jun 2025)
Preprint (discussion started on 17 May 2024)

Interactive discussion

Status: closed

RC1:
'Comment on essd-2024-109', Anonymous Referee #1, 17 Jun 2024

Review comments for the manuscript "Optimal feature selection for improved ML based reconstruction of Global Terrestrial Water Storage Anomalies" of Mandal et al.:
The manuscript describes the derivation of a global long-term terrestrial water storage anomalies data-set, derived from a blending of GRACE satellite observations and global land surface models based on Bayesian networks and machine learning methods.
Such long term information about terrestrial water storage variations is valuable and can help to assess long term trends and to localize extremes. Therefore, I think the manuscript is relevant for publication, however in its current form it lacks to address uncertainties of the product and it is rather structured like a classic scientific study than a data description paper. The structure of the dataset itself is not suitable for efficient usage in the current form and needs revision. Thus for a publication in ESSD, my concerns are as follows.
1) From the title it does not become clear what dataset you would like to advertise. Should it be the TWSAs or the optimal features? I guess its the TWSAs that you would finally like to advertise so you should find a new title in the sense of "ML-based ... long-term terrestrial water storage anomalies from satellite and land-surface model data ...". The term "Optimal feature selection" does not generate an association with a data product, at least for me.
2) In the abstract, you write that you reconstruct TWSA but you don't specify the a grid type and the spatial resolution of the produced dataset. Do you only provide the gridded dataset or also basin aggregates? You should provide this information already in the abstract, although very briefly, so that the reader knows what to expect. Further, I suggest to add a section for data description that explains the structure and content of the final data product in the repository
3) The term "optimal predictors" is mentioned in the title and introduction but the explanation in the methods section (3, 3.1) is not fully clear. What are the optimal predictors? Are they a subsection of your full predictors list? Do you drop training data sets? Are the optimal predictors the ones that have the maximum impact (weight) in the ML algorithms? This should be made more clear and the benefit of knowing the optimal predictors should be outlined.
4) You are not always consistent with your vocabulary, in 3.1 you introduce the term features for what you named previously predictors. I think you should keep a single notion here (and mention the term feature only once, maybe in brackets if this is needed because it's well known by the community).
5) Reproducibility: for making the creation of your data set reproducible, you should at least mention which software tools you used for the machine learning and eventually publish the configurations alongside with your datatset or in another DOI based repository.
6) The selection of evaluation metrics may not be ideal for the global evaluations. CC will always be high for regions with a clear annual amplitude whereas for the deserts with less variations and seasonality it is hard to get a good score in CC. NSE is especially designed for assessing peak flows. Maybe KGE would suit better here. And wouldn't a directed error metric like ME provide additional insight on over- or underestimation tendencies?
7) I think the introduction and the 4.6. Section could be shortened a bit in favor of a data(set) description section.
8) For many of the references DOIs are missing. For several DOI links are incorrect with duplicates in their URLs.
9) Your results should be evaluated in the light of uncertainties of GRACE based water storage anomalies, e.g., https://agupubs.onlinelibrary.wiley.com/doi/10.1029/2021JB022081; there are several different GRACE solutions available which have different levels of uncertainty (https://doi.org/10.1029/2023JB026908) so why did you select exactly one of them and how would the uncertainties of the GRACE product propagate into your BN_TWSA product? The characterization of uncertainties of your gridded and aggregated data sets would be important with respect to deriving any long term trends.
10) The structure of the published dataset does not follow any data standards. Further the naming of the downloadable zip file, Mississipi_Data.zip does not comply with the contained global grids. You should use descriptive filenames and use modern standard data formats, e.g., CDF conform (netCDF) self-describing data, geotiff, ... and a self describing tree structure. You can get inspiration for instance from other publications in ESSD

Minor things:
L86/87: no commas in large numbers
Table 1: Provide not only publications for the data-sets but also the DOI references where they can be obtained; add the acronyms / abbreviations that you later use in the analysis and figures (e.g. Fig. 2, NTWSA, CTWSA)
L199: you name three types of ML algorithms but then are 4 listed and described
L274: That's the third different usage of P in the manuscript (Probability, Precipitation, and Prediction in the evaluation metrics)
L279: A grid is usually defined as a collection of adjacent pixels. You are using the term grid instead of pixel. I suggest to change it to either pixel or grid cell / cells.
Fig.2 Expand acronyms in the figure caption, make the caption more explanatory. From the colors it appears that several optimal predictors overlap for the same regions / pixels
Fig.3 Avoid red and green in the same figure (colorblind check, you can use https://www.color-blindness.com/coblis-color-blindness-simulator/ for checking)
L333: grid-based -> pixel based
Fig. 6: describe gray bars in figure caption (gaps in GRACE solutions)
Fig. 7: Change 1:1 line to non-dashed gray with thicker linewidth to make it distinguishable from the data. Use colors with better contrast for BNML and CTSWA
Abstract L18: remove "and updates will be published when needed"

Citation: https://doi.org/10.5194/essd-2024-109-RC1
- AC1: 'Comment on essd-2024-109', Kironmala Chanda, 08 Jan 2025
  
  The comment was uploaded in the form of a supplement: https://essd.copernicus.org/preprints/essd-2024-109/essd-2024-109-AC1-supplement.pdf
  
  Citation: https://doi.org/10.5194/essd-2024-109-AC1
RC2:
'Comment on essd-2024-109', Anonymous Referee #2, 03 Sep 2024

General comments

Mandal et al. describes a data-driven reconstructed global product of TWS, namely BNML_TWSA. A Bayesian Network technique is used to find an optimal set of predictors. Then, several ML algorithms are trained at each grid cell. The authors choose the best ML algorithm for each grid cell to fore- and hind-cast TWS across the globe. Compared with several existing reconstructed TWS products and estimates by land surface models, the new product shows better agreements with GRACE observations at grid cell and basin scales, with being capable to capture some historical hydroclimatic extreme events. Based on the evaluation results, the authors conclude that the newly developed TWS product is reliable and can be used for hydroclimatological studies.

Overall, I find the manuscript easy to follow. The topic is relevant to the journal, as TWS is an influential variable in many aspects of the Earth system functioning. The evaluation across space and time done by the authors is informative to see how good the ML-derived product performs. However, I still have some concerns regarding the robustness of the new product, majorly coming from the lack of (source of) uncertainty and the way it is evaluated. Please find the comments below. I hope they are helpful to improve the manuscript.

Major comments:

- The the (potential) source(s) of uncertainty needs to be discussed, which the current version of manuscript lacks. There could be some potential sources of uncertainty. First, I wonder if there is any overfit issue in the product as it is fully based on ML and the learning has been done for each grid cell in the study domain. In addition, I wonder if there is any way for BNML_TWSA to capture human impact on TWS and its trend. If no, then it can either be regarded as a source of uncertainty and be specified or be precluded from the model training. If yes, then the authors may add explanations. Lastly, as BNML_TWSA highly depends on TWSA estimates by selected LSMs, the common errors by LSMs, such as the phase shift in mean seasonal cycle (e.g., Bibi et al., 2024) or worse performance in (semi-) arid regions, could propagate to the results. The performance of BNML_TWSA can also potentially be influenced by the precipitation product used.

Bibi, S., Zhu, T., Rateb, A., Scanlon, B. R., Kamran, M. A., Elnashar, A., Bennour, A., and Li, C.: Benchmarking multimodel terrestrial water storage seasonal cycle against Gravity Recovery and Climate Experiment (GRACE) observations over major global river basins, Hydrol. Earth Syst. Sci., 28, 1725–1750, https://doi.org/10.5194/hess-28-1725-2024, 2024.

- I think that readers or potential data users can benefit from additional evaluations. It can be seen that the better performance of BNML_TWSA is expected, as 1) it uses TWSA from LSM(s) as a predictor for many grid cells, 2) the authors choose the best ML algorithms among trained and tested for each grid, and 3) the model is evaluated using TWSA time series. For example, the GRACE-REC by Humphrey and Gudmundsson (2019) that the authors used in the comparison is calibrated against the detrended and deseasonalized TWSA time series, which therefore may not be as good as BNML_TWSA by the design. The results from the evaluation can be seen as the strength of BNML_TWSA, but BNML_TWSA can also benefit from being fairly evaluated against variables that it is not trained with. Evaluating BNML_TWSA using independent variables is important to prove it's ability to extrapolate, because the product has already learned partly the GRACE TWSA information via CTWSA which assimilates GRACE TWSA observations (this fact should have been noted in the main text, I think). As a possible, but not nessesarily the only, way for the additional evaluation, one can suggest that the authors can repeat the evaluation done by Humphrey and Gudmundsson (2019). In the paper, GRACE-REC is evaluated with several independent datasets including sea level budget, streamflow measurements, and basin-scale water balance. This repetition can also work well to compare BNML_TWSA with GRACE-REC. Another way can be evaluating BNML_TWSA at seasonal and interannual (i.e., detrended and deseasonalized) temporal scales. What makes the evaluation at these temporal scales good is because the temporal scales are of a strong interest in multiple communities (e.g., carbon cycle, hydrology, and climate), so the results can be informative for both the product itself and potential users; also it can be regarded as a more fair way to examine BNML_TWSA as it is not trained at this temporal scales.

- I think that the title is misleading. What is the role of the feature selection (i.e., BN) on BNML_TWSA? The title can be seen that using the optimal set of feature given by BN is the key to improve the ML based product in the study, but the relevant section or explanation cannot be seen from the current version of manuscript. So, the contribution of deploying the BN technique to the quality of BNML_TWSA could be more elaborated, or the title could be updated.

- What is the value of examing multiple algorithms for each grid cell? It's clearly reported in the manuscript (e.g., Figure 3) that the spatial pattern leader model is very heterogeneous. However, it has not been reported how different the performance of tested models are, and what the differences are in the actual estimates. I wonder the actual influence on the resulted time series would be minor (e.g., a comparison between the current BNML_TWSA and another BNML_TWSA using the algorithm with the poorest performance for each grid cell), as many ML models usually show similarly good performance.

- The results and discussion section is mostly about presenting how good the performance metrics for BNML_TWSA is, which could have been deeper to provide more insight about the product's applicability. Having evaluation from more diverse aspect as in the second bullet point would help with improving this aspect.

- The dataset provided includes estimates for the Mississippi river basin only, while one would expect estimates for the whole global land grid cells, according to the title and abstract.

Minor and technical comments:

- The current manuscript has many in-text narrative citations that are wrongly used, e.g., L311, L334-335, and so on.

- The introduction can be improved by better addressing the motivation to have a new product. Currently, it introduces TWS variable and examples of (ML-based) TWS reconstruction studies. The authors introduce that using a feature selection process can be a novelty of BNML_TWSAm, while testing multiple ML algorithms is important, which can be, but are not necessarily the reason to have a new product. The authors could better present the motivation by showing why having (or lacking) feature selection and multiple ML algorithms are critical for users and their science. Or, showing from which aspect existing reconstructed TWS products are less reliable/robust can better show the motivation. This will also help with having a focused presentation in the results section.

Specific comments

Abstract

Introduction

L23: There can be more references, especially ones done at the global scale.

L24: This sentence needs more clarification. Which physical processes are missing? What are the influence on the estimates from which aspect?

L29: I think that Mo et al. (2022) is not a proper reference for the sentence. The study is to report a new product, not to examine the human and climatic impact on water cycle.

L41: The reconstruction by Humphrey et al. (2017) is at the global scale. Only the example application is for the Amazon Basin.

L32-59: This paragraph is basically list up previous studies wtih a few sentences for each. I wonder if this is the best way of storytelling for readers.

L60: The authors could first list up what the categories are.

L77-78: As mentioned above, there need to be more elaboration on the importance of the feature selection on the TWSA reconstruction and the applications.

Data and Processing

L108: This counters the sentence L78-80

L112-113: Are there any reasons to choose Noah and CLSM specifically?

L125-126: Why doesn't LSMs fully use the information? This can be more elaborated.

L129: I think that the time period of analysis hasn't specified before.

L143: Although CLSM provides TWS directly, the authors should be able to refer to other materials to know which processes CLSM accounts for to calculate TWS. Please add this description, also plase update Table 1 accordingly.

Eq.1: Sun et al (2019) mentioned that Noah does not account for surface water storage. Please clarify this.

L149: "may be" or "can be"? It's a bit weird to use "may be" in this sentence.

L152-154: It was not clear for me if the prior months are for P and T or for TWSA, too. I suggest to rephrase the sentence.

L153: Why aren't the prior months used for the climate indices?

L155: SVR and ETR need to be introduced as their full name.

Figure 1: For the correlation coefficient, CC has been used through out the manuscript, instead of R.

L231: 'built', not 'build', I think

L233: Is the feature selection procudure in ETR independent to BN?

L266-267: For each grid cell, does BN give the optimal set of predictors to each ML algorithm or is the set common for all ML algorithms? If it's the former, how can one be sure that different ML algorithms share the same optimal set of predictors?

Eq.5: On the right hand side, does the denominator use Pi-Pbar or Oi-Obar? The typical NSE equation uses observations for the denominator. Please check this.

Results and Discussions

Figure 2: Would it be reasonable to interpret the results as the emerging importance of the variables to the global TWSA? For example, can the results be seen as that North Atlantic Oscillation has the least influence on the collective gridwise TWSA among the three modes of climate variability?

L291: Please see the comment for Figure 1

Figure 3: Is it expected that ETR would pop up as the leader? Is there any possible explanation for this, based on the nature of each algorithm?

L306-307: It is not clear that BNML_TWSA performs better than LSMs, especially for the case of CTWSA. Could be clearer with histogram of metrics or the map of differences in metrics between BNML_TWSA and LSMs.

Figure 4: It should be mentioned that BNML_TWSA also shows significant biases in cold or arid regions, and even in a wet region (e.g., a part of the Congo Basin).

Figure 6: What does the shaded area stand for? What is the rationale that BNML_TWSA captures the GRACE TWSA trend, for example, in Indus and GBM Basins, where the TWSA trend would largely be affected by human activity?

Section 4.5: It is great to prove the ability of BNML_TWSA to capture the hydrological extreme events that the MLs were not informed with. However, this section is only mentioning a specific type of event, flood. One could compare the historical TWSA time series of several basins in different climate zones with the corresponding time series of climate indices, drought indices, or precipitation.

Figure 9: I feel that the way to show the ability of BNML_TWSA to capture the historical flood events can be seen as inappropriate. For example, in the map of the USA, I can see many other grid cells as bluish as ones in the green box. Does it mean that all the grid cells similarly bluish as ones in the green box are flooded?

Section 4.6: It is recommended to compare BNML_TWSA and previous studies using independent data sets (please see the first major comment). Also, it needs to be noted that GRACE-REC is calibrated against detrended and deseasonalized GRACE TWSA.

(may not be good at captureing the local TWSA correctly, if GLDAS is calibrated using TWS? BNML_TWSA largely depends on GLDAS LSMs which cannot be reliable at finer spatial scales than the original GRACE spatial resolution) <-- need to check how GLDAS LSMs simulates TWS

Conclusions

L432-434: I think that this point can be mentioned in the main text, possibly with a deeper discussion (i.e., implications of the distribution of selected gridwise predictors for the glbola TWSA and a possible explanation).

Data

- Please provide the unit and the file naming convention.

- This should be also noted in the main text with the number or portion of the grid cells: "# For very few grids BN is not identifying any predictors or only one predictor. Grids with zero or one optimal predictor identified by BN are trained using all 15 potential predictors ('P', 'T', 'NOAH_TWSA', 'CLSM_TWSA', 'DMI', 'NAO', 'ONI', 'P1', 'T1', 'NOAH_TWSA_1', 'CLSM_TWSA_1', 'P2', 'T2', 'NOAH_TWSA_2', 'CLSM_TWSA_2')." Also, I would exclude these grid cells with zero or one predictor identified from Figure 2, if it makes a significant difference.

Citation: https://doi.org/10.5194/essd-2024-109-RC2
- AC1: 'Comment on essd-2024-109', Kironmala Chanda, 08 Jan 2025
  
  The comment was uploaded in the form of a supplement: https://essd.copernicus.org/preprints/essd-2024-109/essd-2024-109-AC1-supplement.pdf
  
  Citation: https://doi.org/10.5194/essd-2024-109-AC1
AC1: 'Comment on essd-2024-109', Kironmala Chanda, 08 Jan 2025

The comment was uploaded in the form of a supplement: https://essd.copernicus.org/preprints/essd-2024-109/essd-2024-109-AC1-supplement.pdf

Citation: https://doi.org/10.5194/essd-2024-109-AC1

Peer review completion

AR: Author's response | RR: Referee report | ED: Editor decision | EF: Editorial file upload

AR by Kironmala Chanda on behalf of the Authors (08 Jan 2025) Author's response Author's tracked changes Manuscript

ED: Referee Nomination & Report Request started (11 Jan 2025) by Kirsten Elger

RR by Anonymous Referee #1 (23 Jan 2025)

Suggestions for revision or reasons for rejection

The authors have satisfactorily revised the manuscript, improved readability, and clarified the raised questions. I have a few more minor remarks that should be addressed before final publication.

1) The findings of the uncertainty assessment (Sec. 4.8) should be shortly reflected in the conclusions as well as in the abstract.

2) Fig. 3: Do you only consider CCs > 0? Otherwise you may find positive differences also by deteriorating from a negative value toward zero.

3) It seems very strange to me that your KGE values are always negative while you get positive NSE values in the range of 0.5 to 0.8. I suggest to reassess the computation of KGE or to provide some discussion why these negative values would be feasible. Possibly, KGE' can be used as an alternative to KGE. Please provide the reference for KGE.

Additional editorial remarks:

L7: I suggest to change grid cell-based to "raster-based"

Table 3 & 4: Captions should not be capitalized

Figure 4: Units are missing for the axes.

Figure 8: Y-Axis labels should be put on the left side, next to the numbers.

Figure 9, 11, 19: Units (mm) are missing for the y-axes.

Figure 18: Notation (Label) for CC is missing for the color-scale

Data product:

Once you publish the global data set, I encourage you to provide also the global grids as netcdf files – instead of only catchment based portions. But probably this is your intention anyway. I think it would be good to provide also the Model_Input and BN_Predictors data as netCDF grids. For the BN_Predictors, you could use an integer based classification. With your current method, you will have to generate a huge number of text-files for the global grid. Importing these data into models will be much easier based on netCDF. Please check also the CF compliance of your netCDF files. Lat and lon seem to have issues here so that cdo does not recognize the grid type (latlong).

Hide

RR by Anonymous Referee #2 (31 Jan 2025)

Suggestions for revision or reasons for rejection

General comments

I appreciate efforts by the authors to address my comments substantially. Now I only have some minor comments.

Minor and technical comments:

- L37-46: I suggest to add a topic sentence at the head of the paragraph, like in the next paragraph. Currently, the paragraph is just a list of example studies and the reader should guess what the paragraph is about based on them.
- L66-87 (A) and L88-109 (B): Each of two paragraphs deals with two topics, so it is less clear what each paragraph wants to talk about. Paragraph A deals with the importance of 1) testing multiple ML algoriths and 2) the optimal feature selection. Paragraph B deals with the importance of the optimal feature selection and debriefs the whole study. I would suggest merge the parts about the optimal feature selection as a new paragraph so that there are three paragraphs about multiple ML algorithms, the optimal feature selection, and debriefing the study.
- Equation 1: I think Eq. 1 can be missled in the current form and a reader may regard that Noah considers surface water storage. Sun et al. (2019) used the equation to show the TWS calculation in general, not for a specific product or model. However, in this study, Equation 1 is to introduce how TWS is calculated in Noah, which does not calculate SWS. So, I suggest to replace CSWS with CWS.
- Figure 4: The six grid cells highlighted are all neighboring to each other. I understand the logic to select the grid cells based on the differences in CC, but using the resulting sample group is less informative for the purpose of the figure, I would say. I suggest to consider selecting other grid cells from other regions, too.
- L358-360: Two last sentences are duplicated in terms of the message.
- L417-418: I think it's wrong to say like "the proposed model accurately reproduced the TWSA series during the pre-GRACE period" by comparing three model simulations.
- L435-436: I suggest to specify the period.
- L515-516: It's hard to agree. Both BNML_TWSA and GRACE struggle to get the observation in rather drier basins, and they behave similarly in general in other basins. One can at least highlight the comparable behavior of Q_BNML_TWSA, or there would need to be additional information (e.g., metrics) to say that Q_BNML_TWSA performs better clearly.
- L521-524: Is it the uncertainty of GRACE as an assumed ground truth that would afffect the product evaluation, or is the processing errors in GRACE liked to the uncertainty of BNML_TWSA? If the former, it would be better to specify it. If it is the latter, BNML_TWSA is still subject to uncertainties of mascon solutions, although they perform better compared to the spherical harmonic solutions. I would indicate both aspects.
- L526-527: The improvement by multiple MLs is majorly for (semi-)arid regions (Figure 3). This means that, every ML algorithm perform greatly for wetter regions or none of them could not improve the simulation there. So, for wetter regions, I suspect there are still rooms to improve (e.g., the KGE map in Figure 6).
- L556-557: The sigma in arid regions may be small by nature because of the small variance of TWS there. However, I agree that testing multiple ML algorithms partly contributed to the improvement in performance (and small sigma values) in arid regions, referring to Figure 3.
- Please make sure that the uncertainty data from Figure 20b is included when the published datasets are updated.

Hide

ED: Publish subject to minor revisions (review by editor) (04 Feb 2025) by Kirsten Elger

AR by Kironmala Chanda on behalf of the Authors (14 Feb 2025) Author's response Author's tracked changes Manuscript

ED: Publish subject to minor revisions (review by editor) (18 Feb 2025) by Kirsten Elger

AR by Kironmala Chanda on behalf of the Authors (21 Feb 2025) Author's response Author's tracked changes Manuscript

ED: Publish as is (04 Mar 2025) by Kirsten Elger

AR by Kironmala Chanda on behalf of the Authors (14 Mar 2025)

Post-review adjustments

AA: Author's adjustment | EA: Editor approval

AA by Kironmala Chanda on behalf of the Authors (02 Jun 2025) Author's adjustment Manuscript

EA: Adjustments approved (03 Jun 2025) by Kirsten Elger

Short summary

Optimal features among hydroclimatic variables and land surface model (LSM) outputs are selected using a novel Bayesian network (BN) approach for simulating terrestrial water storage anomalies (TWSAs). TWSAs are reconstructed (BNML_TWSA) with grid-specific leader models (among four machine learning models) from January 1960 to December 2022 to generate a continuous global gridded dataset. The uncertainty in the reconstructed BNML_TWSA product is also assessed in terms of standard error.