Topology enables learning-based hydrodynamic prediction of the global river system
Authors: Hancheng Ren, Gang Zhao, Shuo Wang, Louise Slater, Dai Yamazaki, Shu Liu, Jingfang Fan, Xueying Li, Shibo Cui, Ziming Yu, Shengyu Kang, Depeng Zuo, Dingzhi Peng, Zongxue Xu, Bo Pang
Organizations: College of Water Sciences, Beijing Normal University, Beijing, China. · School of Geography and the Environment, University of Oxford, Oxford, UK. · Department of Transdisciplinary Science and Engineering, Institute of Science Tokyo, Tokyo, Japan. · School of Systems Science, Beijing Normal University, Beijing, China. · Institute of Industrial Science, University of Tokyo, Tokyo, Japan. · China Institute of Water Resources and Hydropower Research, Beijing, China. · State Key Laboratory of Hydro-Science and Engineering, Department of Hydraulic Engineering, Tsinghua University, Beijing, China. · School of Artificial Intelligence, Beijing Normal University, Beijing, China. · School of Water Resources and Hydropower Engineering, Wuhan University, Wuhan, China.
Accurate river prediction is essential for water, food and energy security, yet remains challenging across entire river networks. Machine learning has transformed Earth-system modeling, but a system-level advance for river prediction lags for lack of reliable data. Exploiting the connectivity and dissipative dynamics of rivers, we introduce GraphRiverCast, a neural model for global river systems that predicts daily multivariate hydrodynamics at every reach of a 0.25°network with only sparse gauges and no initial state. It shows no intrinsic skill decay with lead time, outperforms leading global river models by 25% in accuracy, and robustly generalizes to ungauged reaches and finer resolutions. GraphRiverCast lifts learning-based river prediction from isolated basins to a unified global system and offers insights for machine learning in data-scarce Earth systems.
Figures & tables
Figure 1: Design principle, architecture and training strategy of GraphRiverCast. ( A ) The global river network as a unified system. Each reach carries geomorphic features, runoff forcing, topology and dynamic states, namely discharge ( Q ), water depth ( H ) and storage ( S ). At ungauged reaches (green), no historical states exist, so the standard autoregressive paradigm cannot be initialized. ( B ) Dynamical rationale and architecture. Left, the influence of the initial condition x0 decays as trajectories converge toward a common forcing-driven trajectory after spin-up. Reach-scale prediction can thus be treated as a boundary-value problem, which makes the initial-state-free ColdStart regime physically possible. Center, GraphRiverCast couples a feature encoder, a graph encoder (spatial recursion) and a temporal encoder (temporal recursion) into a topology-informed river model. Right, the fused representation predicts Q , H and S at every reach. ( C ) Training strategy. GraphRiverCast is trained in two stages. It is pre-trained on physical-model simulations over the full network topology, then fine-tuned on sparse in situ stations (red). The fine-tuned model is applied by zero-shot inference across the global river network, shown as a single-time-step discharge snapshot.
Figure 2: Feasibility of initial-state-free prediction and attribution analysis. ( A ) Median combined NSE (mean Nash-Sutcliffe efficiency across discharge, water depth and storage), scored against the CaMa-Flood reference simulation, versus lead time for GRC-ColdStart (blue), GRC-HotStart (red) and a topology-free LSTM run initial-state-free (gray). Bands, 25th–75th percentiles across reaches. Main panel, to a 30-day lead, where the medians cross at day 13 and HotStart falls to 0.40 by day 30 while ColdStart holds at 0.93. Inset, medians over the full test period to a 700-day lead, with the main-panel range shaded and the axis clipped at zero. ( B ) Encoder dependence as barycentric coordinates, the sensitivity of each regime to suppressing each encoder at inference normalized to the share of skill it carries. ( C ) Information flow from inputs to encoders attributed by integrated gradients [ 40 ] , with width proportional to contribution. Only HotStart receives a river-state input. ( D ) Long-lead (21–30 d) NSE stratified by upstream catchment area ( log10 km 2 ) and Shreve magnitude, medians with interquartile bands. At the largest rivers HotStart falls to ≈ 0.30, ColdStart holds ≈ 0.88.
Figure 3: GraphRiverCast benchmarked against leading physics-based and learning-based global river models. ( A ) Global distribution of per-gauge Nash-Sutcliffe efficiency (NSE) for GraphRiverCast over the held-out evaluation period (2018–2019), at the n=1,996 evaluated in situ gauges. Marker color denotes NSE (teal for higher skill, brown for negative), over the global river network (light blue). ( B to E ) Per-gauge comparison of four metrics across the same gauges for GraphRiverCast, an LSTM benchmark, and the physics-based GloFAS and CaMa-Flood, shown as raincloud plots, each combining jittered points (one per gauge), a box for the interquartile range (IQR) with a median line, whiskers at 1.5× IQR and a half-violin kernel density. The value above each distribution is its median, and in E both the median and IQR of PBIAS are given. Higher values indicate better performance for NSE and KGE, lower values for NRMSE, and values near zero for PBIAS. An arrow in each panel marks this direction, upward in B for NSE and C for the Kling-Gupta efficiency (KGE), downward in D for the min-max normalized RMSE, and toward zero in E for the percent bias (PBIAS).
Figure 4: Spatial and scale generalization of GraphRiverCast. ( A ) Per-gauge skill on ungauged reaches, from five-fold cross-validation over the 1,996 evaluable GRDC gauges. The four boxes compare the LSTM and GRC, each trained from scratch or pre-trained then fine-tuned, across four metrics (NSE, KGE, NRMSE, PBIAS). Brackets give differences in medians against the pre-trained LSTM ( Δ Topology), GRC from scratch ( Δ Pre-training) and the LSTM from scratch ( Δ Combined). Asterisks give the significance of the paired differences (Wilcoxon signed-rank test, n.s. not significant). ( B ) The LamaH river network over the Danube at 15 arc-min (214 reaches, top) and 6 arc-min ( 1,409 reaches, bottom). Reaches are blue where resolved at 15 arc-min and orange where resolved only at 6 arc-min. Gauges are red where they map onto a coarse-network reach and gray where recovered only at 6 arc-min. ( C ) Median NSE over the n=325 gauged reaches of the finer network versus the percentage of gauges used for fine-tuning, for GRC fine-tuned from the global pre-trained weights (green), GRC from scratch (gray, dashed) and the pre-trained LSTM (orange). Bands give the min–max across three subsampling seeds, and the filled regions mark where each factor dominates.
Max Planck Institute for Biogeochemistry, Germany · European Centre for Medium-Range Weather Forecasts (ECMWF), United Kingdom · Institute of Bio- and Geosciences Agrosphere (IBG-3), Research Centre Jülich, Germany +2
William & Mary, Williamsburg, VA, USA · Virginia Institute of Marine Science, William & Mary, Gloucester Point, VA, USA · Florida Atlantic University, Boca Raton, FL, USA