the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 04 Aug 2026
RiverGraphNet: Physics-Aware Routing of Gridded Runoff Through Directed River Networks
Abstract. River routing provides the critical link between runoff generation and downstream streamflow prediction, yet conventional routing models often rely on simplified hydraulic assumptions, fixed parameters, and conservative transport formulations that may limit performance in heterogeneous river systems. Here, we introduce RiverGraphNet, a physics-aware graph-based routing framework designed to route physically generated gridded runoff through directed river networks. The framework explicitly isolates routing from runoff generation by coupling Noah-MP runoff with a directed river-network graph derived from the NextGen hydrofabric for the Salt–Verde watershed (Arizona, USA). Gridded runoff is transferred from the Noah-MP domain to graph nodes representing hydrologic routing elements, while streamflow propagation is learned using a Graph Attention Network (GAT) informed by physically meaningful node and edge attributes describing drainage structure, terrain, and hydraulic-routing proxies.
RiverGraphNet was evaluated against observed daily streamflow at 23 USGS gauges and compared with RAPID, a widely used physics-based routing benchmark forced with identical Noah-MP gridded runoff inputs. Experiments examined the influence of temporal routing memory (3–30 day runoff lags), temporal convolution, and alternative loss functions (MSE, weighted MSE, and JKGE-based objectives). RiverGraphNet consistently outperformed RAPID across nearly all gauges and configurations. The best-performing experiment (14-day lag weighted MSE) achieved a median KGE{ss} of 0.72, substantially exceeding RAPID performance (median KGE{ss} of 0.0). Event-scale hydrographs and flow-duration analyses demonstrated improved representation of peak timing, event magnitude, and long-term streamflow distributions. Attention analysis further revealed that routing improvements were most strongly associated with channel-width and dominant-pathway metrics rather than travel-time proxies alone, suggesting that adaptive representation of network influence provides predictive value beyond conventional travel-time parameterization. These results demonstrate the potential of physics-aware, topology-constrained graph learning as a flexible alternative for routing gridded hydrologic runoff through complex river networks.
- Preprint
(6466 KB) - Metadata XML
-
Supplement
(2111 KB) - BibTeX
- EndNote
Status: final response (author comments only)
-
RC1: 'Comment on egusphere-2026-4157', Anonymous Referee #1, 08 Aug 2026
-
AC1: 'Reply on RC1', Mohammad Farmani, 30 Sep 2026
We sincerely thank the referee for the careful evaluation of our manuscript and for the constructive and detailed comments. We appreciate the time and effort devoted to reviewing our work. The comments have helped us improve the clarity of the manuscript, strengthen the methodological discussion, and better define the scope and limitations of the study.
We have carefully considered all comments and revised the manuscript accordingly. In particular, we have clarified the comparison between RiverGraphNet and RAPID, refined the interpretation of model generalization and physical information, added additional analyses addressing multicollinearity and training-seed variability, and expanded the discussion of methodological limitations and reproducibility.
Below, we provide a point-by-point response to each comment. The referee’s comments are reproduced first, followed by our responses and a description of the corresponding revisions made in the manuscript.
Major comments
- The RAPID vs. RiverGraphNet comparison is not on equal footing.
RiverGraphNet is trained in a fully supervised manner directly against observed discharge at the gauges (Section 2.8). RAPID, as described, is not locally recalibrated for this basin, but instead relies on a pre-existing CONUS-scale Muskingum parameterization (as the manuscript itself acknowledges in Section 3.5, lines 753–758). This means the comparison does not isolate "routing representation" alone, as repeatedly stated (e.g., lines 518–521), but conflates two distinct effects: (a) the routing architecture and (b) local calibration against the observed gauge record. The reported RAPID performance (mean KGEss = −0.29, minimum −4.88) is low enough to suggest a parameterization problem rather than an intrinsic limitation of the Muskingum approach. I recommend that the authors: (i) discuss this asymmetry explicitly as a central methodological caveat rather than a passing remark, and (ii) if feasible, include a version of RAPID with at least minimal local calibration (e.g., of the Muskingum X and K parameters) as a fairer benchmark or a sensitivity check.
Response:
We thank the reviewer for raising this important point. We agree that RiverGraphNet and RAPID are not calibrated on exactly equal footing. RiverGraphNet is trained directly against observed streamflow within the Salt–Verde watershed, whereas the RAPID simulation used in this study was not recalibrated specifically for this basin. We have revised the manuscript to make this asymmetry a more explicit methodological caveat and have moderated statements implying that differences between the two models isolate routing architecture alone.
At the same time, we would like to clarify that the RAPID parameterization used here is not an entirely uncalibrated or observation-independent configuration. The adopted parameters originate from a CONUS-scale RAPID calibration framework in which Muskingum parameters were optimized using observed USGS streamflow at multiple gauges. Thus, although we did not perform an additional Salt–Verde-specific recalibration for the present experiment, the regional parameterization was developed using streamflow information from the broader CONUS gauge network rather than purely theoretical or default Muskingum parameters. Large-scale RAPID calibration studies explicitly estimate routing parameters against USGS/NWIS discharge observations, with calibration performed using multiple stream gauges within the river network.
We also agree that the relatively low RAPID KGESS values should not be interpreted as evidence of an intrinsic limitation of the Muskingum method. We have therefore revised the manuscript accordingly. We additionally evaluated RiverGraphNet against the National Water Model (NWM), which uses an extensive calibration and regionalization procedure based on observed streamflow. For the CONUS implementation, NWM parameters are directly calibrated at more than 1,300 relatively unregulated basins using observed discharge, and optimized parameters are subsequently regionalized throughout the remaining domain (Cosgrove et al., 2024). Despite this calibration, NWM performance in the Salt–Verde evaluation was generally lower than that of the RAPID configuration used here, which is why RAPID was retained as the primary and more challenging physics-based benchmark. As shown in Figure 6, RiverGraphNet outperformed both RAPID and NWM across most evaluation gauges.
We therefore interpret the comparison more cautiously in the revised manuscript. The results demonstrate that the supervised RiverGraphNet framework provides substantial predictive improvement relative to the existing large-scale RAPID parameterization and to the calibrated NWM benchmark for this watershed. They do not establish that the GNN architecture is intrinsically superior to a locally optimized Muskingum model. A locally recalibrated RAPID experiment would provide a useful additional benchmark for separating the effects of model architecture from basin-specific parameter optimization, and we now identify this explicitly as an important direction for future evaluation.
To further assess the sensitivity of the comparison to systematic discharge-volume mismatch, we added a training-derived multiplicative correction to RAPID discharge. Correction factors were estimated separately for each gauge using training observations and then held fixed during testing. This increased RAPID mean KGEss from −0.289 to 0.286 and median KGEss from −0.012 to 0.450. The selected RiverGraphNet configuration retained a mean KGEss of 0.696 and a median of 0.723, outperforming corrected RAPID at 21 of 23 gauges. Thus, a constant discharge-volume adjustment reduced the performance gap but did not eliminate RiverGraphNet’s predictive advantage. This experiment adjusts RAPID output magnitude; it does not recalibrate the Muskingum K and X parameters or replace the locally calibrated routing experiment suggested by the reviewer. The procedure and results are reported in Section 2.9 and Section 3.1, respectively, with the comparison summarized in Table 3.
We have clarified the benchmark calibration in Section 2.3 and the evaluation scope in Section 2.9, and discuss the implications of calibration differences in Sections 4.2 and 4.6.
Reference:
Cosgrove, B., Gochis, D., Flowers, T., Dugger, A., Ogden, F., Graziano, T., Clark, E., Cabell, R., Casiday, N., Cui, Z., Eicher, K., Fall, G., Feng, X., Fitzgerald, K., Frazier, N., George, C., Gibbs, R., Hernandez, L., Johnson, D., Jones, R., Karsten, L., Kefelegn, H., Kitzmiller, D., Lee, H., Liu, Y., Mashriqui, H., Mattern, D., McCluskey, A., McCreight, J. L., McDaniel, R., Midekisa, A., Newman, A., Pan, L., Pham, C., RafieeiNasab, A., Rasmussen, R., Read, L., Rezaeianzadeh, M., Salas, F., Sang, D., Sampson, K., Schneider, T., Shi, Q., Sood, G., Wood, A., Wu, W., Yates, D., Yu, W., and Zhang, Y.: NOAA's National Water Model: Advancing operational hydrology through continental‐scale modeling, JAWRA Journal of the American Water Resources Association, 60, 10.1111/1752-1688.13184, 2024.
- Train/test split and generalization claims.
It is not clear whether all 23 gauges appear in both training and test sets (split only by non-overlapping temporal blocks), or whether some gauges are held out entirely to test spatial generalization. If, as it appears, the model observes every gauge-node during training, the strong performance demonstrates temporal interpolation at known nodes rather than generalization to unseen nodes or sub-basins. This should be stated explicitly, since it materially affects the interpretation of claims about "topology-aware learning" and constrains the scope of the transferability discussion in Section 4.6 (point 2).
Response:
We thank the reviewer for identifying this ambiguity. The reviewer is correct that the train, validation, and test partitions in the present study are separated temporally rather than spatially. The same set of gauge locations is represented across the temporal partitions, while non-overlapping 180-day blocks are assigned to training, validation, and testing. Therefore, the reported test performance evaluates generalization to unseen temporal periods at known gauge locations and does not constitute a test of spatial generalization to unseen gauges or sub-basins.
We have revised Section 2.8 to state this explicitly and have clarified the interpretation of the evaluation throughout the manuscript. In particular, we removed language that could imply demonstrated spatial transferability and now distinguish between performance across spatially heterogeneous locations within the study network and transfer to previously unseen network locations.
We have also expanded Section 4.6 to identify spatial generalization as an important limitation of the current study. Future evaluation will require spatial holdout experiments, such as leave-one-gauge-out or leave-subbasin-out testing, as well as multi-basin transfer-learning experiments, to determine whether the topology-aware architecture can generalize to ungauged nodes or river networks not represented during training.
- Use of the term "physics-aware".
The title and abstract emphasize the "physics-aware" nature of the framework, yet physical information (channel width, roughness, etc.) enters only as features of a learned attention mechanism, without any explicit physical constraint (neither mass conservation nor monotonicity with respect to hydraulically meaningful variables). I suggest softening this terminology — e.g., "physically-informed" or "topology-aware with physical attributes" — reserving "physics-aware" for frameworks with explicit constraints, consistent with the authors' own discussion in Sections 4.4/4.6.
Response:
We thank the reviewer for this important clarification. We agree that the original term “physics-aware” could imply a stronger degree of physical constraint than is implemented in the present RiverGraphNet framework. RiverGraphNet incorporates physical information through directed river-network topology and hydraulic, geomorphic, and physiographic node and edge attributes, but it does not currently impose explicit governing equations, mass-conservation constraints, or monotonic hydraulic relationships during learning.
We have therefore revised the terminology throughout the manuscript, replacing “physics-aware” with “physically informed” and, where appropriate, “topology-aware.” The manuscript title has been changed to “RiverGraphNet: Physically Informed Graph-Based Routing of Gridded Runoff Through Directed River Networks.” We have also added explicit clarification in the Methods that “physically informed” refers to the incorporation of physical network structure and attributes rather than explicit enforcement of physical conservation laws.
This revision is consistent with the discussion in Sections 4.4 and 4.6, where we explicitly acknowledge that the present framework does not enforce mass conservation and identify conservation-constrained or hybrid physics–machine-learning routing as an important direction for future work.
- Multicollinearity among attention-related features.
The correlations reported in Table 4 and discussed in Sections 3.7/4.5 between attention coefficients and channel width, travel-time proxy, etc. are interesting, but these features (channel width, travel-time proxy, conveyance proxy, storage proxy) are likely correlated with one another. Without a collinearity analysis (e.g., variance inflation factors or a correlation matrix among the predictors themselves), it is difficult to conclude with confidence that channel width specifically — rather than a correlated feature — is driving attention. I recommend adding this analysis, or at least discussing it explicitly as a limitation.
Response:
We thank the reviewer for this important point. We agree that several of the hydraulic attributes examined in the attention analysis are physically related and may be substantially correlated, which limits attribution of the observed attention relationships to any single predictor.
To address this concern, we performed a collinearity analysis among the principal physical routing predictors used in the attention interpretation. The new Spearman correlation matrix shows substantial dependence among several variables. For example, channel top width is strongly correlated with the storage proxy (r = 0.82), the travel-time proxy is also strongly correlated with the storage proxy (r = 0.82), reach length is correlated with the storage proxy (r = 0.73), and channel slope is negatively correlated with the travel-time proxy (r = -0.71). The maximum variance inflation factor was approximately 9.92, further indicating meaningful multicollinearity among some predictors.
We have added this collinearity analysis to the Supporting Information as Figure S6 and revised Sections 3.7 and 4.5 accordingly. In the revised manuscript, we no longer interpret the relatively strong association between attention-weighted channel-width metrics and routing-performance improvement as evidence that channel width independently drives the learned attention behavior. Instead, we interpret channel width as one indicator within a correlated set of channel-geometry, storage, and hydraulic-conveyance characteristics associated with dominant routing pathways.
We also added this issue explicitly to Section 4.6 as a limitation of the current attention interpretation. The revised manuscript now states that the present analysis identifies statistical associations but cannot uniquely resolve the independent contribution of individual physical attributes. We note that controlled feature-ablation experiments, conditional permutation importance, or other multivariate attribution methods would be required to isolate these effects more rigorously.
- Statistical robustness of the GNN results.
It is not stated whether the results in Table 3 and Figure 5 derive from a single training run or from multiple runs with different random seeds. Neural networks carry non-negligible stochastic variability from initialization and training; without a confidence interval on the model's KGEss (analogous to the standard deviation reported across gauges for RAPID), the Wilcoxon test (lines 556–559) captures variability across gauges but not variability arising from the model's own training process.
Response:
We thank the reviewer for highlighting this issue. Each RiverGraphNet configuration was trained using multiple random seeds, and the original manuscript reported the best-performing realization for each configuration. We agree that this selection procedure was not sufficiently documented and that the cross-gauge variability reported in Table S2 (formerly Table 3) does not quantify stochastic uncertainty associated with neural-network initialization and optimization.
We have therefore added a multi-seed robustness analysis in the revised manuscript. For each training realization, we calculated the mean test across evaluation gauges and summarized its variability across random seeds. The new Figure S2 and Table S3 report the distribution of seed-level performance, mean, standard deviation, and range across successfully evaluated seeds.
The analysis reveals an additional result that was not apparent from the selected best-seed experiments. Several CNN-enabled configurations exhibit substantial sensitivity to random initialization, with some seeds converging to high-performing solutions and others producing near-zero or negative . In contrast, successful no-CNN configurations show markedly lower seed-to-seed variability. For example, for the lag14 weighted-MSE experiments, the CNN-enabled configuration achieved a mean gauge-averaged of 0.47 with an across-seed standard deviation of 0.35, whereas the no-CNN configuration achieved 0.67 with a standard deviation of 0.01 among valid runs. Similar behavior was observed for the MSE and JKGE objectives. We have revised Section 3.3 to discuss this additional evidence that the temporal convolutional head is not necessary for strong routing skill and may increase optimization sensitivity.
Table S2 reports the selected best-performing realization for each configuration, whereas the new supplementary analysis quantifies stochastic training variability across seeds. We have also clarified that the paired Wilcoxon test assesses consistency of performance differences across gauges for the selected model realization and does not quantify training-seed uncertainty; the latter is now assessed separately through the multi-seed analysis.
- The "JKGE" metric and its citation.
The "Jawad Kling–Gupta Efficiency" metric is named after a co-author of the present manuscript (Muhammad Jawad) and cites a work that appears to be unpublished or in press in the same year. Please verify that this citation is complete and verifiable by reviewers and readers, and consider whether the naming convention could give an impression of insufficiently transparent self-citation.
Response:
We thank the reviewer for raising this concern. The manuscript cited the JKGE formulation incompletely. The associated study is currently under review at Water Resources Research, but a publicly accessible preprint is available on arXiv (Jawad et al., 2026, arXiv:2604.03906). We have revised the reference to provide the complete and verifiable arXiv citation.
We also agree that repeatedly referring to the metric as the “Jawad Kling–Gupta Efficiency” could create unnecessary ambiguity because Muhammad Jawad is also a co-author of the present manuscript. We have therefore revised the terminology throughout the manuscript to refer to the metric primarily as JKGE, with attribution to Jawad et al. (2026), and have added a brief methodological description explaining that it is a modification of KGE designed to account for temporal non-stationarity in the benchmark statistical properties.
These revisions make clear that the metric was developed in a separate, publicly available study and is not introduced or named within the present RiverGraphNet manuscript.
Minor comments
- Line 194: "(Noah-MP, citation?)" is an unresolved placeholder and must be fixed before publication.
Response: Fixed
- Line 599: typo "occurrs" → "occurs".
Response: Fixed
- Lines 707/711 and nearby: several subject–verb agreement slips ("reproduces" where "reproduce" is needed); a full language edit is recommended.
Response: Fixed
- Figure 8: the legend includes a ">300%" improvement class, while the text (lines 512–514) explicitly cautions that percent improvements should be interpreted with care when the RAPID denominator is near zero. Consider removing this class from the figure or flagging in the caption which points fall into this problematic regime.
Response: We thank the reviewer for highlighting this issue. Figure 8 displays percentage improvement relative to RAPID, and we agree that very large percentage values can arise when the corresponding RAPID is close to zero. We retained the percentage representation because it provides a useful spatial visualization of relative performance gains but revised the figure caption and Section 3.2 to explicitly caution that the > 300% category may reflect sensitivity to a small RAPID denominator and should not be interpreted as a precise measure of improvement magnitude. The figure is therefore used primarily to identify the spatial distribution of relative gains rather than to emphasize the exact magnitude of extreme percentage values.
- Computational cost (training time, parameter count) of the GNN framework is not reported; this would be a useful practical comparison against RAPID, which is far cheaper to run (though not to calibrate).
Response:
We thank the reviewer for this suggestion. The selected RiverGraphNet configuration contained 17,842 trainable parameters, which we have now reported in Section 2.7. In our implementation, training required approximately 8 h on a GPU, with variability across configurations and random seeds, while inference required less than 10 min. These timings describe our implementation rather than a controlled runtime comparison with RAPID. RAPID was executed on CPUs and RiverGraphNet on GPUs, and their workflows include different preprocessing and calibration requirements. We therefore report the model size and approximate computational requirements without drawing a direct conclusion about their relative computational efficiency.
- It is unclear whether the Noah-MP/AORC forcing data and RAPID outputs are publicly available alongside the GitHub code repository; this would improve full reproducibility.
Response: We thank the reviewer for pointing out this ambiguity. We have clarified the Data and Code Availability section to distinguish between publicly available source datasets and the model-generated products used in this study. The AORC meteorological forcing is publicly available, whereas the Noah-MP runoff simulations and RAPID routing simulations were performed offline. Because the complete forcing, Noah-MP runoff, and RAPID output datasets are very large, they are not included in the GitHub repository. The repository instead provides the RiverGraphNet source code, training configurations, graph-cache files, evaluation scripts, and supporting analysis workflows used in the graph-routing experiments.
- Section 2.7: Equations 1–2 are correct, but the number of GAT layers, hidden dimension, and dropout used in the main experiment should be stated in the text itself, not only in Figure 1.
Response: We thank the reviewer for this suggestion. We have revised Section 2.7 to explicitly state the principal GAT architectural settings used in the main experiments. The RiverGraphNet architecture consisted of four GAT layers with a hidden feature dimension of 64, four attention heads, and a dropout rate of 0.1. These settings were previously shown only in Figure 1 and are now also reported in the Methods text for clarity and reproducibility.
Major comments
- The RAPID vs. RiverGraphNet comparison is not on equal footing.
RiverGraphNet is trained in a fully supervised manner directly against observed discharge at the gauges (Section 2.8). RAPID, as described, is not locally recalibrated for this basin, but instead relies on a pre-existing CONUS-scale Muskingum parameterization (as the manuscript itself acknowledges in Section 3.5, lines 753–758). This means the comparison does not isolate "routing representation" alone, as repeatedly stated (e.g., lines 518–521), but conflates two distinct effects: (a) the routing architecture and (b) local calibration against the observed gauge record. The reported RAPID performance (mean KGEss = −0.29, minimum −4.88) is low enough to suggest a parameterization problem rather than an intrinsic limitation of the Muskingum approach. I recommend that the authors: (i) discuss this asymmetry explicitly as a central methodological caveat rather than a passing remark, and (ii) if feasible, include a version of RAPID with at least minimal local calibration (e.g., of the Muskingum X and K parameters) as a fairer benchmark or a sensitivity check.
Response:
We thank the reviewer for raising this important point. We agree that RiverGraphNet and RAPID are not calibrated on exactly equal footing. RiverGraphNet is trained directly against observed streamflow within the Salt–Verde watershed, whereas the RAPID simulation used in this study was not recalibrated specifically for this basin. We have revised the manuscript to make this asymmetry a more explicit methodological caveat and have moderated statements implying that differences between the two models isolate routing architecture alone.
At the same time, we would like to clarify that the RAPID parameterization used here is not an entirely uncalibrated or observation-independent configuration. The adopted parameters originate from a CONUS-scale RAPID calibration framework in which Muskingum parameters were optimized using observed USGS streamflow at multiple gauges. Thus, although we did not perform an additional Salt–Verde-specific recalibration for the present experiment, the regional parameterization was developed using streamflow information from the broader CONUS gauge network rather than purely theoretical or default Muskingum parameters. Large-scale RAPID calibration studies explicitly estimate routing parameters against USGS/NWIS discharge observations, with calibration performed using multiple stream gauges within the river network.
We also agree that the relatively low RAPID KGESS values should not be interpreted as evidence of an intrinsic limitation of the Muskingum method. We have therefore revised the manuscript accordingly. We additionally evaluated RiverGraphNet against the National Water Model (NWM), which uses an extensive calibration and regionalization procedure based on observed streamflow. For the CONUS implementation, NWM parameters are directly calibrated at more than 1,300 relatively unregulated basins using observed discharge, and optimized parameters are subsequently regionalized throughout the remaining domain (Cosgrove et al., 2024). Despite this calibration, NWM performance in the Salt–Verde evaluation was generally lower than that of the RAPID configuration used here, which is why RAPID was retained as the primary and more challenging physics-based benchmark. As shown in Figure 6, RiverGraphNet outperformed both RAPID and NWM across most evaluation gauges.
We therefore interpret the comparison more cautiously in the revised manuscript. The results demonstrate that the supervised RiverGraphNet framework provides substantial predictive improvement relative to the existing large-scale RAPID parameterization and to the calibrated NWM benchmark for this watershed. They do not establish that the GNN architecture is intrinsically superior to a locally optimized Muskingum model. A locally recalibrated RAPID experiment would provide a useful additional benchmark for separating the effects of model architecture from basin-specific parameter optimization, and we now identify this explicitly as an important direction for future evaluation.
To further assess the sensitivity of the comparison to systematic discharge-volume mismatch, we added a training-derived multiplicative correction to RAPID discharge. Correction factors were estimated separately for each gauge using training observations and then held fixed during testing. This increased RAPID mean KGEss from −0.289 to 0.286 and median KGEss from −0.012 to 0.450. The selected RiverGraphNet configuration retained a mean KGEss of 0.696 and a median of 0.723, outperforming corrected RAPID at 21 of 23 gauges. Thus, a constant discharge-volume adjustment reduced the performance gap but did not eliminate RiverGraphNet’s predictive advantage. This experiment adjusts RAPID output magnitude; it does not recalibrate the Muskingum K and X parameters or replace the locally calibrated routing experiment suggested by the reviewer. The procedure and results are reported in Section 2.9 and Section 3.1, respectively, with the comparison summarized in Table 3.
We have clarified the benchmark calibration in Section 2.3 and the evaluation scope in Section 2.9, and discuss the implications of calibration differences in Sections 4.2 and 4.6.
Reference:
Cosgrove, B., Gochis, D., Flowers, T., Dugger, A., Ogden, F., Graziano, T., Clark, E., Cabell, R., Casiday, N., Cui, Z., Eicher, K., Fall, G., Feng, X., Fitzgerald, K., Frazier, N., George, C., Gibbs, R., Hernandez, L., Johnson, D., Jones, R., Karsten, L., Kefelegn, H., Kitzmiller, D., Lee, H., Liu, Y., Mashriqui, H., Mattern, D., McCluskey, A., McCreight, J. L., McDaniel, R., Midekisa, A., Newman, A., Pan, L., Pham, C., RafieeiNasab, A., Rasmussen, R., Read, L., Rezaeianzadeh, M., Salas, F., Sang, D., Sampson, K., Schneider, T., Shi, Q., Sood, G., Wood, A., Wu, W., Yates, D., Yu, W., and Zhang, Y.: NOAA's National Water Model: Advancing operational hydrology through continental‐scale modeling, JAWRA Journal of the American Water Resources Association, 60, 10.1111/1752-1688.13184, 2024.
- Train/test split and generalization claims.
It is not clear whether all 23 gauges appear in both training and test sets (split only by non-overlapping temporal blocks), or whether some gauges are held out entirely to test spatial generalization. If, as it appears, the model observes every gauge-node during training, the strong performance demonstrates temporal interpolation at known nodes rather than generalization to unseen nodes or sub-basins. This should be stated explicitly, since it materially affects the interpretation of claims about "topology-aware learning" and constrains the scope of the transferability discussion in Section 4.6 (point 2).
Response:
We thank the reviewer for identifying this ambiguity. The reviewer is correct that the train, validation, and test partitions in the present study are separated temporally rather than spatially. The same set of gauge locations is represented across the temporal partitions, while non-overlapping 180-day blocks are assigned to training, validation, and testing. Therefore, the reported test performance evaluates generalization to unseen temporal periods at known gauge locations and does not constitute a test of spatial generalization to unseen gauges or sub-basins.
We have revised Section 2.8 to state this explicitly and have clarified the interpretation of the evaluation throughout the manuscript. In particular, we removed language that could imply demonstrated spatial transferability and now distinguish between performance across spatially heterogeneous locations within the study network and transfer to previously unseen network locations.
We have also expanded Section 4.6 to identify spatial generalization as an important limitation of the current study. Future evaluation will require spatial holdout experiments, such as leave-one-gauge-out or leave-subbasin-out testing, as well as multi-basin transfer-learning experiments, to determine whether the topology-aware architecture can generalize to ungauged nodes or river networks not represented during training.
- Use of the term "physics-aware".
The title and abstract emphasize the "physics-aware" nature of the framework, yet physical information (channel width, roughness, etc.) enters only as features of a learned attention mechanism, without any explicit physical constraint (neither mass conservation nor monotonicity with respect to hydraulically meaningful variables). I suggest softening this terminology — e.g., "physically-informed" or "topology-aware with physical attributes" — reserving "physics-aware" for frameworks with explicit constraints, consistent with the authors' own discussion in Sections 4.4/4.6.
Response:
We thank the reviewer for this important clarification. We agree that the original term “physics-aware” could imply a stronger degree of physical constraint than is implemented in the present RiverGraphNet framework. RiverGraphNet incorporates physical information through directed river-network topology and hydraulic, geomorphic, and physiographic node and edge attributes, but it does not currently impose explicit governing equations, mass-conservation constraints, or monotonic hydraulic relationships during learning.
We have therefore revised the terminology throughout the manuscript, replacing “physics-aware” with “physically informed” and, where appropriate, “topology-aware.” The manuscript title has been changed to “RiverGraphNet: Physically Informed Graph-Based Routing of Gridded Runoff Through Directed River Networks.” We have also added explicit clarification in the Methods that “physically informed” refers to the incorporation of physical network structure and attributes rather than explicit enforcement of physical conservation laws.
This revision is consistent with the discussion in Sections 4.4 and 4.6, where we explicitly acknowledge that the present framework does not enforce mass conservation and identify conservation-constrained or hybrid physics–machine-learning routing as an important direction for future work.
- Multicollinearity among attention-related features.
The correlations reported in Table 4 and discussed in Sections 3.7/4.5 between attention coefficients and channel width, travel-time proxy, etc. are interesting, but these features (channel width, travel-time proxy, conveyance proxy, storage proxy) are likely correlated with one another. Without a collinearity analysis (e.g., variance inflation factors or a correlation matrix among the predictors themselves), it is difficult to conclude with confidence that channel width specifically — rather than a correlated feature — is driving attention. I recommend adding this analysis, or at least discussing it explicitly as a limitation.
Response:
We thank the reviewer for this important point. We agree that several of the hydraulic attributes examined in the attention analysis are physically related and may be substantially correlated, which limits attribution of the observed attention relationships to any single predictor.
To address this concern, we performed a collinearity analysis among the principal physical routing predictors used in the attention interpretation. The new Spearman correlation matrix shows substantial dependence among several variables. For example, channel top width is strongly correlated with the storage proxy (r = 0.82), the travel-time proxy is also strongly correlated with the storage proxy (r = 0.82), reach length is correlated with the storage proxy (r = 0.73), and channel slope is negatively correlated with the travel-time proxy (r = -0.71). The maximum variance inflation factor was approximately 9.92, further indicating meaningful multicollinearity among some predictors.
We have added this collinearity analysis to the Supporting Information as Figure S6 and revised Sections 3.7 and 4.5 accordingly. In the revised manuscript, we no longer interpret the relatively strong association between attention-weighted channel-width metrics and routing-performance improvement as evidence that channel width independently drives the learned attention behavior. Instead, we interpret channel width as one indicator within a correlated set of channel-geometry, storage, and hydraulic-conveyance characteristics associated with dominant routing pathways.
We also added this issue explicitly to Section 4.6 as a limitation of the current attention interpretation. The revised manuscript now states that the present analysis identifies statistical associations but cannot uniquely resolve the independent contribution of individual physical attributes. We note that controlled feature-ablation experiments, conditional permutation importance, or other multivariate attribution methods would be required to isolate these effects more rigorously.
- Statistical robustness of the GNN results.
It is not stated whether the results in Table 3 and Figure 5 derive from a single training run or from multiple runs with different random seeds. Neural networks carry non-negligible stochastic variability from initialization and training; without a confidence interval on the model's KGEss (analogous to the standard deviation reported across gauges for RAPID), the Wilcoxon test (lines 556–559) captures variability across gauges but not variability arising from the model's own training process.
Response:
We thank the reviewer for highlighting this issue. Each RiverGraphNet configuration was trained using multiple random seeds, and the original manuscript reported the best-performing realization for each configuration. We agree that this selection procedure was not sufficiently documented and that the cross-gauge variability reported in Table S2 (formerly Table 3) does not quantify stochastic uncertainty associated with neural-network initialization and optimization.
We have therefore added a multi-seed robustness analysis in the revised manuscript. For each training realization, we calculated the mean test across evaluation gauges and summarized its variability across random seeds. The new Figure S2 and Table S3 report the distribution of seed-level performance, mean, standard deviation, and range across successfully evaluated seeds.
The analysis reveals an additional result that was not apparent from the selected best-seed experiments. Several CNN-enabled configurations exhibit substantial sensitivity to random initialization, with some seeds converging to high-performing solutions and others producing near-zero or negative . In contrast, successful no-CNN configurations show markedly lower seed-to-seed variability. For example, for the lag14 weighted-MSE experiments, the CNN-enabled configuration achieved a mean gauge-averaged of 0.47 with an across-seed standard deviation of 0.35, whereas the no-CNN configuration achieved 0.67 with a standard deviation of 0.01 among valid runs. Similar behavior was observed for the MSE and JKGE objectives. We have revised Section 3.3 to discuss this additional evidence that the temporal convolutional head is not necessary for strong routing skill and may increase optimization sensitivity.
Table S2 reports the selected best-performing realization for each configuration, whereas the new supplementary analysis quantifies stochastic training variability across seeds. We have also clarified that the paired Wilcoxon test assesses consistency of performance differences across gauges for the selected model realization and does not quantify training-seed uncertainty; the latter is now assessed separately through the multi-seed analysis.
- The "JKGE" metric and its citation.
The "Jawad Kling–Gupta Efficiency" metric is named after a co-author of the present manuscript (Muhammad Jawad) and cites a work that appears to be unpublished or in press in the same year. Please verify that this citation is complete and verifiable by reviewers and readers, and consider whether the naming convention could give an impression of insufficiently transparent self-citation.
Response:
We thank the reviewer for raising this concern. The manuscript cited the JKGE formulation incompletely. The associated study is currently under review at Water Resources Research, but a publicly accessible preprint is available on arXiv (Jawad et al., 2026, arXiv:2604.03906). We have revised the reference to provide the complete and verifiable arXiv citation.
We also agree that repeatedly referring to the metric as the “Jawad Kling–Gupta Efficiency” could create unnecessary ambiguity because Muhammad Jawad is also a co-author of the present manuscript. We have therefore revised the terminology throughout the manuscript to refer to the metric primarily as JKGE, with attribution to Jawad et al. (2026), and have added a brief methodological description explaining that it is a modification of KGE designed to account for temporal non-stationarity in the benchmark statistical properties.
These revisions make clear that the metric was developed in a separate, publicly available study and is not introduced or named within the present RiverGraphNet manuscript.
Minor comments
- Line 194: "(Noah-MP, citation?)" is an unresolved placeholder and must be fixed before publication.
Response: Fixed
- Line 599: typo "occurrs" → "occurs".
Response: Fixed
- Lines 707/711 and nearby: several subject–verb agreement slips ("reproduces" where "reproduce" is needed); a full language edit is recommended.
Response: Fixed
- Figure 8: the legend includes a ">300%" improvement class, while the text (lines 512–514) explicitly cautions that percent improvements should be interpreted with care when the RAPID denominator is near zero. Consider removing this class from the figure or flagging in the caption which points fall into this problematic regime.
Response: We thank the reviewer for highlighting this issue. Figure 8 displays percentage improvement relative to RAPID, and we agree that very large percentage values can arise when the corresponding RAPID is close to zero. We retained the percentage representation because it provides a useful spatial visualization of relative performance gains but revised the figure caption and Section 3.2 to explicitly caution that the > 300% category may reflect sensitivity to a small RAPID denominator and should not be interpreted as a precise measure of improvement magnitude. The figure is therefore used primarily to identify the spatial distribution of relative gains rather than to emphasize the exact magnitude of extreme percentage values.
- Computational cost (training time, parameter count) of the GNN framework is not reported; this would be a useful practical comparison against RAPID, which is far cheaper to run (though not to calibrate).
Response:
We thank the reviewer for this suggestion. The selected RiverGraphNet configuration contained 17,842 trainable parameters, which we have now reported in Section 2.7. In our implementation, training required approximately 8 h on a GPU, with variability across configurations and random seeds, while inference required less than 10 min. These timings describe our implementation rather than a controlled runtime comparison with RAPID. RAPID was executed on CPUs and RiverGraphNet on GPUs, and their workflows include different preprocessing and calibration requirements. We therefore report the model size and approximate computational requirements without drawing a direct conclusion about their relative computational efficiency.
- It is unclear whether the Noah-MP/AORC forcing data and RAPID outputs are publicly available alongside the GitHub code repository; this would improve full reproducibility.
Response: We thank the reviewer for pointing out this ambiguity. We have clarified the Data and Code Availability section to distinguish between publicly available source datasets and the model-generated products used in this study. The AORC meteorological forcing is publicly available, whereas the Noah-MP runoff simulations and RAPID routing simulations were performed offline. Because the complete forcing, Noah-MP runoff, and RAPID output datasets are very large, they are not included in the GitHub repository. The repository instead provides the RiverGraphNet source code, training configurations, graph-cache files, evaluation scripts, and supporting analysis workflows used in the graph-routing experiments.
- Section 2.7: Equations 1–2 are correct, but the number of GAT layers, hidden dimension, and dropout used in the main experiment should be stated in the text itself, not only in Figure 1.
Response: We thank the reviewer for this suggestion. We have revised Section 2.7 to explicitly state the principal GAT architectural settings used in the main experiments. The RiverGraphNet architecture consisted of four GAT layers with a hidden feature dimension of 64, four attention heads, and a dropout rate of 0.1. These settings were previously shown only in Figure 1 and are now also reported in the Methods text for clarity and reproducibility.
Citation: https://doi.org/10.5194/egusphere-2026-4157-AC1
-
AC1: 'Reply on RC1', Mohammad Farmani, 30 Sep 2026
-
RC2: 'Comment on egusphere-2026-4157', Anonymous Referee #2, 22 Sep 2026
The authors present a graph-based neural network system for routing gridded runoff as an alternative to conventional routing methods (i.e. Muskingum-based RAPID and NWM) as applied in large-scale hydrologic models. The manuscript describes the method as applied to watersheds in Arizona using runoff from a NOAA-MP model forced with historical meteorology and finds that the neural network approach can achieve much better replication of observed hydrographs at gages throughout the region.
The manuscript is generally written well, is well organized, and is thorough in its coverage of the interpretation and implications of the analysis. The overall effort is an interesting one and has the potential to be of considerable interest to the community, especially with regard to the evaluation of different loss functions configurations and their impact on the results as well as the possible connection between learned routing metrics and performance improvement. However, I am concerned the framing and presentation of the work in its current form is possibly misleading and should be addressed. These concerns and some minor comments are noted below.
Major comments
For much of the manuscript, the extent to which physics are incorporated in the "physics-aware" RiverGraphNet is not very clear. It appears the physics involved are primarily the topology of the river network and the gridded runoff inherited from Noah-MP. Importantly, there is no mass balance constraint on the RiverGraphNet routing. Although this is explained in the later sections of the manuscript, I'd suggest clarifying more specifically and earlier in the manuscript what is and is not explicitly included in the neural network so that the reader can evaluate the results in a more appropriate context. Reframing this as "topologically aware" or something similar would be more appropriate than the current "physics-aware" and related terminology used throughout.
Related to this, the manuscript presents the RiverGraphNet method as superior to traditional routing methods (like RAPID) with the assumption that 1) the Noah-MP runoff data and observed gage data are already consistent and 2) the RAPID method provides a poor match to the gage record because it is inherently limited or flawed. Although the authors do address this later in the manuscript, the bulk of the text treats as secondary the possibility that there are biases in the underlying gridded runoff from Noah-MP that a strictly mass conservative routing approach cannot overcome, regardless of whether the Muskingum routing method is well calibrated or even appropriately suited for parts of the watershed (it might or might not be!). If there is a chance that the Noah-MP runoff estimates have biases in volume or timing (and the NWM/RAPID reference metrics suggest this might be the case) while the objective is to fit observed gage records with no mass balance constraint, then one would naturally expect the unconstrained RiverGraphNet method to perform substantially better. A fairer --and potentially more useful -- comparison would be one that removes any underlying biases in the gridded runoff first before performing the training and analysis. Again, this context needs to be addressed much earlier in the manuscript along with a more detailed justification of the value the current analysis provides.
The mass balance evaluation (comparison of upstream runoff to streamflow volume at a gage) is great and I'm glad it was included. I'd suggest potentially including it in the main portion of the manuscript along with a bit more detail on how this was calculated. For example, are the runoff/streamflow volume ratios compared on a total volume per period basis or as an average of ratios calculated on each day? The former seems like the better approach for an aggregate water budget, but it's unclear if that's the case.
Much of the discussion provided around the results and interpretation alludes to some complexity (seasonal, spatial, antecedent hydrologic context, etc) that fixed parameter routing methods cannot inherently capture. I think it would be valuable to be more specific about what these complexities are, what the current evidence for the region tells us about them, and how the RiverGraphNet helps elucidate possible solutions. For example, is the study characterized by reaches with strongly (and potentially seasonally) losing conditions? Do the results from RiverGraphNet corroborate this? Are channel configurations in some areas substantially different after big runoff/snowmelt events compared to after months/years of drought? Is there a signal from the RiverGraphNet results that can point to where this might occur? Ultimately, the findings that would be most impactful and interesting from this work would be related to how the neural network method can help identify where existing traditional routing methods like RAPID can/should be modified and in what manner. The authors mention how some of these ideas may be applied in future work, but it seems much of the information for laying a more substantial foundation for this is already at hand. A bit more hydrologic context provided at the beginning of the manuscript tied to the interpretation of the results would greatly improve the impact of the paper.
Minor comments
Check that all acronyms are defined in full before use in the main body and abstract. For example, GAT is defined in full in the abstract and then used initially in the main body as an acronym before being defined again.
Table 3 - I am not sure if this table provides much value as Figure 6 serves as an effective visual summary of the same information. Perhaps this can be moved to the supplement.
There are numerous minor instances of unrendered or unresolved formatting throughout - I suggest scanning through to check for and fix such occurrences.
Page 21 - Figure 7 & 8 (and related discussion) - An approximate visual comparison of the information in Figures 7 and 8 suggests that the magnitude of KGEss improvement is inversely correlated with the RAPID KGE: if the RAPID KGE was higher, the improvement with the GAT tends to be lower. If this is the case, then there is perhaps less support for interpretating improvement patterns based on landscape location or related characteristics. Some additional discussion and clarification of this in the manuscript would be helpful.
Citation: https://doi.org/10.5194/egusphere-2026-4157-RC2 -
AC2: 'Reply on RC2', Mohammad Farmani, 30 Sep 2026
We sincerely thank the referee for the careful evaluation of our manuscript and for the constructive and detailed comments. We appreciate the time and effort devoted to reviewing our work. The comments have helped us improve the clarity of the manuscript, strengthen the methodological discussion, and better define the scope and limitations of the study.
Below, we provide a point-by-point response to each comment. The referee’s comments are reproduced first, followed by our responses and a description of the corresponding revisions made in the manuscript.
Major comments
- For much of the manuscript, the extent to which physics are incorporated in the "physics-aware" RiverGraphNet is not very clear. It appears the physics involved are primarily the topology of the river network and the gridded runoff inherited from Noah-MP. Importantly, there is no mass balance constraint on the RiverGraphNet routing. Although this is explained in the later sections of the manuscript, I'd suggest clarifying more specifically and earlier in the manuscript what is and is not explicitly included in the neural network so that the reader can evaluate the results in a more appropriate context. Reframing this as "topologically aware" or something similar would be more appropriate than the current "physics-aware" and related terminology used throughout.
Response:
We thank the reviewer for this important comment. We agree that the original terminology could imply a stronger degree of physical constraint than is implemented in RiverGraphNet. This concern was also raised by another reviewer, therfore we have revised the manuscript throughout accordingly. Specifically, we replaced the term “physics-aware” with “physically informed” and revised the title to “RiverGraphNet: Physically Informed Graph-Based Routing of Gridded Runoff Through Directed River Networks.” We also added explicit clarification in Section 2.7 that the framework incorporates physical information through directed river-network topology and hydraulic, geomorphic, and physiographic attributes, but does not impose explicit governing equations, mass-conservation constraints, momentum conservation, or monotonic hydraulic relationships.
In response to the present reviewer’s suggestion that this distinction should be introduced earlier, we have also added a clarification in Section 1.3 stating that “physically informed” refers to the incorporation of physical network structure and attributes rather than explicit enforcement of physical conservation laws. This allows the scope and limitations of the framework to be clear before the reader reaches the Results and Discussion.
- Related to this, the manuscript presents the RiverGraphNet method as superior to traditional routing methods (like RAPID) with the assumption that 1) the Noah-MP runoff data and observed gage data are already consistent and 2) the RAPID method provides a poor match to the gage record because it is inherently limited or flawed. Although the authors do address this later in the manuscript, the bulk of the text treats as secondary the possibility that there are biases in the underlying gridded runoff from Noah-MP that a strictly mass conservative routing approach cannot overcome, regardless of whether the Muskingum routing method is well calibrated or even appropriately suited for parts of the watershed (it might or might not be!). If there is a chance that the Noah-MP runoff estimates have biases in volume or timing (and the NWM/RAPID reference metrics suggest this might be the case) while the objective is to fit observed gage records with no mass balance constraint, then one would naturally expect the unconstrained RiverGraphNet method to perform substantially better. A fairer --and potentially more useful -- comparison would be one that removes any underlying biases in the gridded runoff first before performing the training and analysis. Again, this context needs to be addressed much earlier in the manuscript along with a more detailed justification of the value the current analysis provides.
Response:
We thank the reviewer for highlighting this important distinction. We agree that using identical Noah-MP runoff inputs does not, by itself, isolate differences in routing representation. RiverGraphNet does not enforce mass conservation and can therefore improve agreement with observed discharge partly through compensation for errors in runoff magnitude or timing. A conservative routing framework cannot correct an upstream volume mismatch through redistribution alone. Consequently, the original performance comparison should not be interpreted as demonstrating the intrinsic superiority of RiverGraphNet’s routing formulation.
To investigate the contribution of systematic volume mismatch, we conducted an additional sensitivity experiment using a training-derived multiplicative correction to RAPID discharge. For each gauge, the correction factor was calculated as the ratio of accumulated reference discharge to accumulated RAPID discharge over matching valid training dates. This factor was then held fixed and applied to RAPID discharge during the held-out test blocks. No validation or test observations were used to estimate the factors. All models were evaluated over the same valid test dates using the discharge reference employed in the manuscript.
The correction substantially improved RAPID performance: the mean KGEss increased from −0.289 to 0.286, and median KGEss increased from −0.012 to 0.450. The selected RiverGraphNet configuration achieved a mean KGEss of 0.696 and a median of 0.723, outperforming volume-corrected RAPID at 21 of the 23 evaluation gauges. Thus, systematic volume mismatch contributes materially to the original performance gap, while a constant volume adjustment does not eliminate RiverGraphNet’s predictive advantage under the evaluation used here.
This additional analysis directly tests whether a constant, training-derived adjustment of discharge volume can reproduce the observed predictive gains. The improvement in RAPID performance demonstrates the importance of systematic volume mismatch, while RiverGraphNet’s retained advantage shows that this adjustment alone does not reproduce its predictive performance. We have incorporated both findings into the manuscript and revised the interpretation of the original benchmark comparison accordingly.
We acknowledge that this experiment is a discharge-volume correction sensitivity analysis, rather than the bias-corrected gridded-runoff experiment suggested by the reviewer. The gauge-specific multipliers do not define a spatially consistent correction across nested catchments, and they do not modify event timing or represent time-varying errors. The remaining performance difference therefore cannot be attributed exclusively to improved routing physics; it may also reflect more flexible compensation for forcing errors and unrepresented hydrologic processes.
A more comprehensive evaluation would apply a spatially consistent runoff correction before routing and reassess both frameworks under the corrected forcing, however, this is beyond the scope of the current manuscript. We are also developing a differentiable Muskingum routing implementation that would enable joint optimization of a runoff-correction network and routing parameters against downstream discharge. This framework remains under development and is intended for a separate study; it is not part of the present RiverGraphNet experiments.
- The mass balance evaluation (comparison of upstream runoff to streamflow volume at a gage) is great and I'm glad it was included. I'd suggest potentially including it in the main portion of the manuscript along with a bit more detail on how this was calculated. For example, are the runoff/streamflow volume ratios compared on a total volume per period basis or as an average of ratios calculated on each day? The former seems like the better approach for an aggregate water budget, but it's unclear if that's the case.
Response:
We thank the reviewer for requesting this clarification. The diagnostic is calculated as total streamflow volume divided by total upstream Noah-MP runoff volume over the same retained dates, rather than as an average of daily ratios. Section 2.9 now explains the unit conversions, upstream aggregation, common valid-date selection, and interpretation of the discontinuous evaluation intervals. We have moved the diagnostic into the main manuscript as Figure 12 and discuss it in Section 3.6. Because storage changes across interval boundaries are not explicitly included, we describe it as a runoff–streamflow consistency diagnostic rather than a water-balance closure test.
The following is added to the section 2.9:
“We additionally evaluated runoff–streamflow consistency using the ratio of total streamflow volume to total upstream Noah-MP runoff volume at each gauge. Daily discharge (m³ s⁻¹) was multiplied by 86,400 s to obtain daily streamflow volume. Daily surface and subsurface runoff depths (mm day⁻¹) were summed, multiplied by the mapped contributing grid-cell areas (m²), and divided by 1,000 to obtain runoff volumes. These contributions were accumulated over all graph nodes upstream of and including the gauge node. For each gauge, streamflow and runoff volumes were summed over the same retained dates before division. Thus, the diagnostic represents a ratio of total volumes, rather than an average of daily ratios. Dates were retained only when observations, RiverGraphNet predictions, RAPID discharge, and aggregated runoff were available and finite.
The diagnostic used the temporally held-out test blocks, with 30 days removed from each end of the overall test-date range. The retained dates therefore comprise separated intervals rather than a continuous evaluation period. Because changes in channel storage across interval boundaries are not explicitly included, this ratio is interpreted as a runoff–streamflow consistency diagnostic rather than a strict water-balance closure test. Departures from unity may reflect runoff-input errors, storage changes, unrepresented water gains or losses, and routing behavior.”
- Much of the discussion provided around the results and interpretation alludes to some complexity (seasonal, spatial, antecedent hydrologic context, etc) that fixed parameter routing methods cannot inherently capture. I think it would be valuable to be more specific about what these complexities are, what the current evidence for the region tells us about them, and how the RiverGraphNet helps elucidate possible solutions. For example, is the study characterized by reaches with strongly (and potentially seasonally) losing conditions? Do the results from RiverGraphNet corroborate this? Are channel configurations in some areas substantially different after big runoff/snowmelt events compared to after months/years of drought? Is there a signal from the RiverGraphNet results that can point to where this might occur? Ultimately, the findings that would be most impactful and interesting from this work would be related to how the neural network method can help identify where existing traditional routing methods like RAPID can/should be modified and in what manner. The authors mention how some of these ideas may be applied in future work, but it seems much of the information for laying a more substantial foundation for this is already at hand. A bit more hydrologic context provided at the beginning of the manuscript tied to the interpretation of the results would greatly improve the impact of the paper.
Response:
We thank the reviewer for encouraging a more specific connection between regional hydrology and the interpretation of our results. We have added regional context in Section 2.1, drawing on studies of the Verde River that document spatially variable groundwater contributions and streamflow losses, including contrasting conditions during the February 2011 and June 2007 surveys (Garner and Bills, 2012), and gaining, losing, and dry reaches in the middle Verde watershed (Paretti et al., 2018). The revised text therefore recognizes both groundwater contributions and streamflow losses rather than characterizing the network uniformly as losing.
In Section 4.5, we connect this evidence to potential extensions of the evaluated routing configuration. The training-derived discharge-volume correction improved RAPID mean test KGEss from −0.289 to 0.286, while RiverGraphNet retained a mean of 0.696 and outperformed corrected RAPID at 21 of 23 gauges. This result demonstrates predictive improvements beyond those obtained through a constant gauge-specific volume adjustment and motivates testing more flexible representations of runoff-to-discharge relationships. We discuss explicit groundwater exchange and channel storage, alongside channel evaporation and riparian evapotranspiration, as candidate processes for future experiments retaining Muskingum routing and accounting for water gains, losses, and storage changes.
The present diagnostics do not independently identify individual losing reaches, quantify exchange fluxes, or establish changes in channel geometry following floods or droughts. The regional literature provides physical motivation for the proposed extensions, while RiverGraphNet provides predictive evidence supporting their investigation. We have revised the discussion to distinguish these hypotheses from processes demonstrated by the present experiments.
References:
Garner, B. D. and Bills, D. J.: Spatial and seasonal variability of base flow in the Verde Valley, central Arizona, 2007 and 2011, i-33, 10.3133/sir20125192, 2012.
Paretti, N. V., Brasher, A. M. D., Pearlstein, S. L., Skow, D. M., Gungle, B. W., and Garner, B. D.: Preliminary synthesis and assessment of environmental flows in the middle Verde River watershed, Arizona, 10.3133/sir20175100, 2018.
Minor comments
- Check that all acronyms are defined in full before use in the main body and abstract. For example, GAT is defined in full in the abstract and then used initially in the main body as an acronym before being defined again.
Response: We have checked acronym definitions in the abstract and main text separately and expanded terms at their first occurrence. We have also standardized subsequent usage, including the definitions of GAT, USGS, RAPID, MSE, and clarified the terminology used for the NextGen framework.
- Table 3 - I am not sure if this table provides much value as Figure 6 serves as an effective visual summary of the same information. Perhaps this can be moved to the supplement.
Response: The original Table 3, which summarized performance across model configurations, has been moved to Table S2. The revised Table 3 presents the newly added discharge-volume-correction sensitivity experiment and therefore contains a different comparison.
- There are numerous minor instances of unrendered or unresolved formatting throughout - I suggest scanning through to check for and fix such occurrences.
Response: We have reviewed the manuscript and Supporting Information and corrected unresolved mathematical text, equation numbering, punctuation, spacing, and figure/table cross-references. We have also checked the presentation of mathematical symbols and subscripts.
- Page 21 - Figure 7 & 8 (and related discussion) - An approximate visual comparison of the information in Figures 7 and 8 suggests that the magnitude of KGEss improvement is inversely correlated with the RAPID KGE: if the RAPID KGE was higher, the improvement with the GAT tends to be lower. If this is the case, then there is perhaps less support for interpretating improvement patterns based on landscape location or related characteristics. Some additional discussion and clarification of this in the manuscript would be helpful.
Response: We thank the reviewer for this observation. We examined the relationship between gauge-level improvement and baseline RAPID performance and confirmed the pattern identified by the reviewer. For the selected 14-day-lag weighted-MSE configuration, absolute improvement in KGEss was strongly negatively associated with RAPID KGEss across the 23 evaluation gauges (Spearman ρ = −0.968). Associations between improvement and upstream drainage area or stream order were substantially weaker (ρ = −0.243 and −0.248, respectively).
This relationship partly reflects the definition of improvement as RiverGraphNet KGEss minus RAPID KGEss, together with the comparatively narrow distribution of RiverGraphNet scores. Therefore, it does not independently demonstrate that particular landscape characteristics explain the performance gains.
We have revised Section 3.2 to report these relationships and distinguish the geographic distribution of predictive improvements from their physical interpretation. We have removed statements attributing larger gains to headwater or intermediate-basin processes based on the maps alone, and clarified that smaller gains may reflect stronger baseline RAPID performance. We also revised the figure captions to emphasize interpretation alongside absolute model skill and the sensitivity of percentage improvements to near-zero RAPID scores.
Citation: https://doi.org/10.5194/egusphere-2026-4157-AC2
-
AC2: 'Reply on RC2', Mohammad Farmani, 30 Sep 2026
Viewed
| HTML | XML | Total | Supplement | BibTeX | EndNote | |
|---|---|---|---|---|---|---|
| 189 | 122 | 29 | 340 | 49 | 21 | 20 |
- HTML: 189
- PDF: 122
- XML: 29
- Total: 340
- Supplement: 49
- BibTeX: 21
- EndNote: 20
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
This manuscript presents a graph attention network (GAT) framework for routing Noah-MP-generated runoff through the NGen river-network topology of the Salt–Verde watershed, and compares it against RAPID. The strategy of isolating the routing step from runoff generation is well motivated, and the set of experiments (temporal-lag configurations, loss-function formulations, ablation of the temporal convolutional head, and attention analysis) is thorough and clearly documented. The candid discussion in Sections 4.4 and 4.6, where the authors acknowledge that the framework does not enforce mass conservation and that part of the performance gain may reflect implicit compensation for Noah-MP biases, is a genuine strength and is not common in manuscripts of this kind.
That said, several issues should be addressed before the manuscript is suitable for publication. These are listed below as major and minor comments.
Major comments
RiverGraphNet is trained in a fully supervised manner directly against observed discharge at the gauges (Section 2.8). RAPID, as described, is not locally recalibrated for this basin, but instead relies on a pre-existing CONUS-scale Muskingum parameterization (as the manuscript itself acknowledges in Section 3.5, lines 753–758). This means the comparison does not isolate "routing representation" alone, as repeatedly stated (e.g., lines 518–521), but conflates two distinct effects: (a) the routing architecture and (b) local calibration against the observed gauge record. The reported RAPID performance (mean KGEss = −0.29, minimum −4.88) is low enough to suggest a parameterization problem rather than an intrinsic limitation of the Muskingum approach. I recommend that the authors: (i) discuss this asymmetry explicitly as a central methodological caveat rather than a passing remark, and (ii) if feasible, include a version of RAPID with at least minimal local calibration (e.g., of the Muskingum X and K parameters) as a fairer benchmark or a sensitivity check.
It is not clear whether all 23 gauges appear in both training and test sets (split only by non-overlapping temporal blocks), or whether some gauges are held out entirely to test spatial generalization. If, as it appears, the model observes every gauge-node during training, the strong performance demonstrates temporal interpolation at known nodes rather than generalization to unseen nodes or sub-basins. This should be stated explicitly, since it materially affects the interpretation of claims about "topology-aware learning" and constrains the scope of the transferability discussion in Section 4.6 (point 2).
The title and abstract emphasize the "physics-aware" nature of the framework, yet physical information (channel width, roughness, etc.) enters only as features of a learned attention mechanism, without any explicit physical constraint (neither mass conservation nor monotonicity with respect to hydraulically meaningful variables). I suggest softening this terminology — e.g., "physically-informed" or "topology-aware with physical attributes" — reserving "physics-aware" for frameworks with explicit constraints, consistent with the authors' own discussion in Sections 4.4/4.6.
The correlations reported in Table 4 and discussed in Sections 3.7/4.5 between attention coefficients and channel width, travel-time proxy, etc. are interesting, but these features (channel width, travel-time proxy, conveyance proxy, storage proxy) are likely correlated with one another. Without a collinearity analysis (e.g., variance inflation factors or a correlation matrix among the predictors themselves), it is difficult to conclude with confidence that channel width specifically — rather than a correlated feature — is driving attention. I recommend adding this analysis, or at least discussing it explicitly as a limitation.
It is not stated whether the results in Table 3 and Figure 5 derive from a single training run or from multiple runs with different random seeds. Neural networks carry non-negligible stochastic variability from initialization and training; without a confidence interval on the model's KGEss (analogous to the standard deviation reported across gauges for RAPID), the Wilcoxon test (lines 556–559) captures variability across gauges but not variability arising from the model's own training process.
The "Jawad Kling–Gupta Efficiency" metric is named after a co-author of the present manuscript (Muhammad Jawad) and cites a work that appears to be unpublished or in press in the same year. Please verify that this citation is complete and verifiable by reviewers and readers, and consider whether the naming convention could give an impression of insufficiently transparent self-citation.
Minor comments
Summary for the authors
This is a methodologically interesting and unusually candid piece of work regarding its own limitations, but Major Comment 1 (the asymmetric comparison between a gauge-calibrated model and an uncalibrated benchmark) is central and needs to be addressed directly. Otherwise, the paper's headline message — that the GNN substantially outperforms RAPID — risks being more misleading than informative.