River Flow and Nutrient Load Projection Guide
river-projection(stage 1, below) lives inocean-prep, alongsidebc-correct/bc-diagnostics— install and run it from theocean-preprepo (cli/river_projection.py), with its config atocean-prep/config/river_projection_example.yaml.river-diagnostics(stage 2) stays inocean-post. Confirmed 2026-10-07; this page previously referenced aconfig/river_projection.yamlinocean-postthat was never there, and the real tool and config were found inocean-prepinstead. No project-specific copy of the config exists yet in this checkout — the only one found isocean-prep’s own example.
Overview
This guide describes how to project future river flows and nutrient loadings using the delta-change method: bias-corrected gridded precipitation (P) and evapotranspiration (E) from CMIP6 are used to scale the observed seasonal discharge cycle from the EMORID dataset forward in time. Nutrient loads are then computed by multiplying projected discharge by an observed flow-weighted mean concentration.
EMORID observed Q + nutrients
│
├─── monthly Q climatology ──────────────────────────┐
│ (calibration period) │
│ ▼
ERA5 P, E (reference) Q_future = Q_clim(month) × ΔPE(t)
│ │
└─── P–E climatology ────────────────► ΔPE(t) │
(calibration period, denominator) │
│
BC future P, E (CMIP6) ▼
│ nutrient loads = Q_future × FWMC
└─── future P–E ──────────────────────► (numerator)
Commands
# Stage 1 (ocean-prep): project future flows and nutrient loads
river-projection --config config/river_projection_example.yaml --model MPI-ESM1-2-HR --scenario ssp245
river-projection --config config/river_projection_example.yaml --model MPI-ESM1-2-HR --scenario ssp585
river-projection --config config/river_projection_example.yaml --model MPI-ESM1-2-HR --scenario ssp245 --dryrun
river-projection --config config/river_projection_example.yaml --model MPI-ESM1-2-HR --scenario ssp245 --allow-missing
# Stage 2 (ocean-post): generate diagnostic plots and tables from the projection output
river-diagnostics --config config/river_projection_example.yaml --model MPI-ESM1-2-HR --scenario ssp245
river-diagnostics --config config/river_projection_example.yaml --model MPI-ESM1-2-HR --scenario ssp245 --analyses-dir /data/local/analyses
--model is required by river-projection (either metadata.model in the
config or --model) — it resolves {model} in the forcing and output paths
alongside --scenario.
Config on GitHub: config/river_projection_example.yaml (in ocean-prep)
Input data
| Dataset | Location | Content |
|---|---|---|
| EMORID observed flows | ${EMORID_FOLDER}/to_netcdf/EMORID_1990_2024.nc | daily Q + nutrients |
| ERA5 precipitation | ${ERA5_FOLDER}/era5_tp_{year}.nc | Hourly accumulation (m), variable tp — reference P–E climatology |
| ERA5 evaporation | ${ERA5_FOLDER}/era5_e_{year}.nc | Hourly accumulation (m), variable e, negative for upward flux |
| BC future P/E (separate) | from bc-correct (ocean-prep) | Daily, one file per year, variables pr/evspsbl (kg m⁻² s⁻¹) |
| BC future P−E (composite) | from bc-correct (ocean-prep) | Jointly corrected pe, preferred over separate pr/evspsbl when present — see “The pe composite” below |
EMORID variables available: Q, TotalN, NO3, NH4, DIN, TotalP, PO4, Si, DOC, DIC, SPM. Countries covered: 18 (NE Atlantic and Baltic rim).
Method
Step 1 — Station filtering
Stations are filtered to retain only those with ≥ min_years (default 10) of
valid daily Q. This removes gauges with very short or mostly-missing records
that would produce unreliable climatologies.
Stations in Norway, Finland, Sweden, and Iceland are flagged with
regulated_flag = True in the output but are not excluded. Rivers in
these countries are heavily influenced by snowmelt timing and hydropower
operations, so the P–E relationship does not hold as reliably as for
rainfall-dominated catchments. Projections for flagged stations should be
interpreted with caution.
Step 2 — Associate stations with gridded P–E
Each station is paired with the nearest grid cell in the ERA5 P and E fields (for the climatology) and the bias-corrected future P and E fields (for the projection) using nearest-neighbour selection on (lat, lon).
This is a practical approximation: in the absence of catchment boundary polygons the point extraction captures the local climate signal but not the integrated catchment response. A future improvement would ingest catchment boundary polygons (e.g. HydroSHEDS) and compute proper area-weighted P–E over each drainage basin. The delta-change method is robust to this approximation because it relies only on the relative change in P–E, not its absolute value.
Step 3 — Monthly P–E climatology (calibration period)
The mean monthly P–E at each station is computed over the calibration period (e.g. 1990–2019) from ERA5, which is the observational reference:
PE_ref_clim(m) = mean over t in calibration period of [P_ERA5(t) - E_ERA5(t)]
where month(t) = m
This 12-value climatology forms the denominator of the delta-change factor.
Why ERA5 and not bias-corrected historical? The BC future files are anchored to the ERA5 climatological distribution by construction. Dividing the future P–E by the ERA5 climatology therefore recovers the pure CMIP6 climate-change signal. Bias-corrected historical files would give the same result (ERA5 ≈ BC-historical by construction) and are not needed.
ERA5 units and sign convention: ERA5 tp/e are hourly accumulations
in metres; CMIP6’s bias-corrected pr/evspsbl are daily-mean fluxes in
kg m⁻² s⁻¹. Both ERA5 variables are resampled to daily totals
(resample: daily, summing the 24 hourly values) and then scaled by
0.011574 (= 1000/86400/24) to match that flux convention. ERA5 e is
also stored as a negative accumulation (upward flux = surface loss), so
its scale factor is -0.011574 — the sign flip and the unit conversion
combined into one number.
Step 4 — Observed monthly Q climatology (calibration period)
The mean monthly observed discharge at each station is computed from EMORID over the same calibration period:
Q_clim(m) = mean over t in calibration period of Q_obs(t)
where month(t) = m
This is the baseline seasonal cycle that the delta-change method preserves.
Step 5 — Delta-change factor
For each future time step t the delta-change factor is:
ΔPE(t) = PE_future(t) / PE_hist_clim(month(t))
where PE_future(t) = P_future(t) − E_future(t) from the bias-corrected
future forcing.
Interpretation: ΔPE(t) = 1.2 means the future P–E in that month is 20 % higher than the historical climatological value, so the projected discharge is also 20 % higher.
Handling non-positive P–E: In arid regions or summer months, the historical P–E climatology can be zero or negative (evaporation exceeds precipitation). Division by a non-positive value is physically undefined for this method; those cells are set to ΔPE = 1.0 (no change) and logged as warnings. If many stations trigger this condition, consider extending the calibration period or reviewing the evaporation variable.
Step 6 — Projected river discharge
Q_future(t) = Q_clim(month(t)) × ΔPE(t)
Values are clipped to zero (no negative discharge).
Why delta-change rather than a full rainfall-runoff model? The delta-change method preserves the observed seasonal cycle and inter-annual variability pattern of discharge without requiring catchment area, soil parameters, or routing. It applies only the climate-driven fractional change in water availability (P–E) to the observed flow. The method is widely used for first-order assessments when catchment characteristics are not available.
Step 7 — Flow-weighted mean concentration (FWMC)
For each nutrient species and station:
FWMC = Σ(Q(t) × C(t)) / Σ(Q(t)) over the FWMC period
The FWMC period is configured independently from the calibration period. Use a recent sub-period (e.g. 2010–2023) rather than the full record to reflect current nutrient conditions; older data from the 1990s may overestimate future loads in rivers that have improved under EU nutrient directives such as the Water Framework Directive and Nitrates Directive.
Both Q and C must be finite for a time step to enter the sum. Stations with
no valid Q–C overlap for a given nutrient receive FWMC = NaN and are
excluded from load projections for that nutrient.
Step 8 — Projected nutrient loads
Load_future(t) = Q_future(t) × FWMC [kg s⁻¹ if C in kg m⁻³]
This assumes in-catchment nutrient sources and in-stream processing remain at their current level. This is a conservative baseline; a separate scenario analysis can relax it by applying trend-based or policy-driven adjustments to FWMC.
Configuration
# ocean-prep/config/river_projection_example.yaml
metadata:
scenario: ssp126 # override at run time with --scenario
model: MPI-ESM1-2-HR # override at run time with --model
emorid:
path: "${EMORID_FOLDER}/to_netcdf/EMORID_1990_2024.nc"
station_filter:
min_years: 10 # minimum years of valid daily Q
calibration:
start: "1990-01-01"
end: "2019-12-31"
future:
start: "2015-01-01"
end: "2099-12-31"
fwmc: # period for flow-weighted mean concentrations
start: "2010-01-01" # omit this block to use the full EMORID record
end: "2023-12-31"
nutrients:
- TotalN
- TotalP
- DIN
domain:
lon_min: -20.0
lon_max: 40.0
lat_min: 45.0
lat_max: 75.0
forcing:
# ERA5: denominator of the delta-change factor (P–E climatology).
reference:
pr:
path: "${ERA5_FOLDER}/era5_tp_{year}.nc"
variable: tp
scale_factor: 0.011574 # hourly accumulation (m) → daily kg m⁻² s⁻¹
resample: daily # hourly → daily sum
evspsbl:
path: "${ERA5_FOLDER}/era5_e_{year}.nc"
variable: e
scale_factor: -0.011574 # sign flip (ERA5 e is negative) + unit conversion
resample: daily
# BC future: numerator of the delta-change factor. {model}/{scenario}
# substituted from metadata/--model/--scenario; {method} (regridding
# method, e.g. qdm/bilinear) is left as a glob (*) rather than hardcoded,
# since pr and evspsbl may use different methods.
future:
pr:
path: "${BIAS_CORRECTED_FOLDER}/CMIP6/{model}/{scenario}/meteo/pr_bc_*_{scenario}_{year}.nc"
variable: pr
evspsbl:
path: "${BIAS_CORRECTED_FOLDER}/CMIP6/{model}/{scenario}/meteo/evspsbl_bc_*_{scenario}_{year}.nc"
variable: evspsbl
# Jointly bias-corrected P-E composite, preferred over separate pr −
# evspsbl when present (auto-detected by real file presence per
# model/scenario — not every model produces it, e.g. GFDL-ESM4's CMIP6
# hfls is Amon-only). See "The pe composite" below.
pe:
path: "${BIAS_CORRECTED_FOLDER}/CMIP6/{model}/{scenario}/meteo/pe_bc_*_{scenario}_{year}.nc"
variable: pe
output:
dir: "${BIAS_CORRECTED_FOLDER}/CMIP6/{model}/{scenario}/rivers"
compress: true
analyses_dir: "${OCEANICU_ANALYSES_FOLDER}" # not used by river-projection itself today —
# a placeholder for river-diagnostics (stage 2)
{model} and {scenario} in the future forcing paths and output.dir are
substituted at runtime from metadata.model/metadata.scenario (or
--model/--scenario). This makes it straightforward to run multiple
models and SSPs from the one shared config:
for scenario in ssp126 ssp370 ssp585; do
river-projection --config config/river_projection_example.yaml \
--model MPI-ESM1-2-HR --scenario $scenario
done
The pe composite
When forcing.future.pe is present in the config, river-projection probes
for real output at that path for the current --model/--scenario before
using it: if found, it uses the jointly bias-corrected P−E composite instead
of separately-corrected pr − evspsbl (separate correction does not
guarantee the difference preserves ERA5’s own P–E climatology). If no
pe output exists for that model/scenario, it logs this and falls back to
pr/evspsbl automatically — the same shared config works for every model
via --model without a per-model variant.
Output
Two NetCDF files are written per run, both with dimensions (site, time):
| File | Variables | Description |
|---|---|---|
river_flows_future_<scenario>.nc | Q, regulated_flag | Projected monthly discharge (m³ s⁻¹) |
nutrient_loads_future_<scenario>.nc | TotalN, TotalP, DIN, … | Projected loads (kg s⁻¹) |
Global attributes record the scenario name, calibration and future periods,
and the FWMC period. These are written to output.dir (not into the
analyses_dir/validations tree) — output.analyses_dir in the config is
not used by river-projection itself; it only matters to river-diagnostics
below.
Diagnostics
After running river-projection, use river-diagnostics (ocean-post) to
generate plots and tables from the projection output. It reads the same
config file (for metadata, emorid, and output.dir) and writes into
the canonical analyses layout:
<analyses_dir>/scenarios/<model>/<scenario>/river/
plots/
q_change_map_<scenario>.png # station map coloured by % Q change
q_time_series_<scenario>.png # annual mean Q: calibration + near/mid/far future
q_seasonal_<scenario>.png # monthly Q climatology by period
nutrient_loads_<scenario>.png # annual nutrient load time series
tables/
river_projection_statistics_<scenario>.csv
river_projection_statistics_<scenario>.txt
river-diagnostics --config config/river_projection_example.yaml --model MPI-ESM1-2-HR --scenario ssp245
river-diagnostics --config config/river_projection_example.yaml --model MPI-ESM1-2-HR --scenario ssp585
river-diagnostics --config config/river_projection_example.yaml --model MPI-ESM1-2-HR --scenario ssp245 --analyses-dir /data/local/analyses
--analyses-dir overrides output.analyses_dir from the config and has no
./analyses fallback — see
Tidal Analysis: analyses_dir has no default.
--scenario and --model override the corresponding metadata keys.
The diagnostics are automatically picked up by ocean-reporting and included
in the Hugo scenarios page under a River Flow and Nutrient Load Projections
section.
Dry run
Use --dryrun to check the config and probe forcing-file existence (every
probed year, for both ERA5 reference and BC future — plus the pe auto-
detect probe when configured) without loading any data:
river-projection --config config/river_projection_example.yaml --model MPI-ESM1-2-HR --scenario ssp245 --dryrun
Both --dryrun and a real run exit non-zero if any probed year is missing —
a config/filename mismatch is meant to fail loudly, not silently degrade the
projection. Pass --allow-missing only for genuinely partial, expected
coverage (e.g. a scenario whose later years are not produced yet).
Limitations and future improvements
| Limitation | Impact | Potential improvement |
|---|---|---|
| Nearest-neighbour P–E extraction | Ignores catchment shape; grid noise can affect individual stations | Ingest HydroSHEDS catchment polygons and compute area-weighted P–E |
| Constant FWMC | Does not capture land-use change or emission trends | Fit a trend to historical C(t) and extrapolate; or apply emission scenario adjustments |
| Regulated rivers flagged but not excluded | Delta-change projections less physically meaningful for snowmelt / hydropower dominated basins | Separate treatment using temperature-based snowmelt model |
| No uncertainty from FWMC | Single FWMC value per station | Bootstrap confidence intervals on FWMC from the observed record |
Key references
- Arnell (1999), Climate change and global water resources, Global Environmental Change 9, S31–S49. — Foundation of the delta-change approach for river flow projections.
- Middelkoop et al. (2001), Impact of climate change on hydrological regimes and water resources management in the Rhine basin, Climatic Change 49, 105–128. — Applied delta-change example.
- Kronvang et al. (2009), Nitrogen and phosphorus losses from agricultural areas in European river basins, Science of the Total Environment. — FWMC methodology context for riverine nutrient loads.