Working with Ensembles¶

Every generator and disaggregator in SynHydro returns an Ensemble: a container that holds many synthetic realizations for many sites, plus metadata about how they were produced. This tutorial covers the data model, how to build an Ensemble from your own data, and the built-in analysis and I/O methods.

Note: Run this notebook from the examples/ directory. Figures are saved to figures/04_ensembles/ and data outputs to outputs/04_ensembles/ (both are gitignored).

In [1]:

from pathlib import Path

FIG_DIR = Path("figures") / "04_ensembles"
OUT_DIR = Path("outputs") / "04_ensembles"
FIG_DIR.mkdir(parents=True, exist_ok=True)
OUT_DIR.mkdir(parents=True, exist_ok=True)

from pathlib import Path FIG_DIR = Path("figures") / "04_ensembles" OUT_DIR = Path("outputs") / "04_ensembles" FIG_DIR.mkdir(parents=True, exist_ok=True) OUT_DIR.mkdir(parents=True, exist_ok=True)

The dual data model¶

An Ensemble maintains two complementary views of the same data:

View	Keys	Each value
data_by_realization	realization index (int)	DataFrame [time x sites]
data_by_site	site name (str)	DataFrame [time x realizations]

Use whichever fits the task: per-realization DataFrames feed simulation models one trace at a time, while per-site DataFrames make it easy to compute distributional statistics across realizations.

In [2]:

import synhydro

Q_daily = synhydro.load_example_data()
Q_monthly = Q_daily.resample("MS").sum()

gen = synhydro.KirschGenerator()
gen.fit(Q_monthly)
ensemble = gen.generate(n_realizations=50, n_years=30, seed=42)

print(ensemble.realization_ids[:5], "...")     # [0, 1, ..., 49]
print(ensemble.site_names)                     # site names from the fitted data

trace = ensemble.data_by_realization[0]                # one trace, all sites
site_flows = ensemble.data_by_site["USGS-01434000"]    # one site, all realizations
trace.head()

import synhydro Q_daily = synhydro.load_example_data() Q_monthly = Q_daily.resample("MS").sum() gen = synhydro.KirschGenerator() gen.fit(Q_monthly) ensemble = gen.generate(n_realizations=50, n_years=30, seed=42) print(ensemble.realization_ids[:5], "...") # [0, 1, ..., 49] print(ensemble.site_names) # site names from the fitted data trace = ensemble.data_by_realization[0] # one trace, all sites site_flows = ensemble.data_by_site["USGS-01434000"] # one site, all realizations trace.head()

[0, 1, 2, 3, 4] ...
['USGS-01434000', 'USGS-01438500', 'USGS-01440000', 'USGS-01463500']

Out[2]:

	USGS-01434000	USGS-01438500	USGS-01440000	USGS-01463500
2026-01-01	2195.882255	2400.959069	50.496287	5716.950581
2026-02-01	6585.782036	7207.159678	127.695394	13044.662260
2026-03-01	5388.598558	6523.424728	110.861873	12214.155599
2026-04-01	6538.286359	7317.730322	147.966110	15561.684423
2026-05-01	2939.916869	3808.781743	96.352207	10024.805864

The by-site view is computed lazily from the by-realization view the first time you access it, so holding an Ensemble costs no more memory than the realization dictionary until you need the second view.

Building an Ensemble from your own data¶

Pass a dictionary in either orientation; the structure is detected from the key type. Integer keys mean the dictionary is keyed by realization, string keys mean it is keyed by site. First build some small DataFrames to use as stand-ins for your own data:

In [3]:

import numpy as np
import pandas as pd

rng = np.random.default_rng(1)
years = pd.date_range("2000-01-01", periods=10, freq="YS")

# Two realizations, each [time x sites]
df_realization_0 = pd.DataFrame(
    rng.gamma(2.0, 100.0, size=(10, 2)), index=years, columns=["Hopland", "Healdsburg"]
)
df_realization_1 = pd.DataFrame(
    rng.gamma(2.0, 100.0, size=(10, 2)), index=years, columns=["Hopland", "Healdsburg"]
)

# Two sites, each [time x realizations]
df_hopland = pd.DataFrame(rng.gamma(2.0, 100.0, size=(10, 3)), index=years, columns=[0, 1, 2])
df_healdsburg = pd.DataFrame(rng.gamma(2.0, 100.0, size=(10, 3)), index=years, columns=[0, 1, 2])

import numpy as np import pandas as pd rng = np.random.default_rng(1) years = pd.date_range("2000-01-01", periods=10, freq="YS") # Two realizations, each [time x sites] df_realization_0 = pd.DataFrame( rng.gamma(2.0, 100.0, size=(10, 2)), index=years, columns=["Hopland", "Healdsburg"] ) df_realization_1 = pd.DataFrame( rng.gamma(2.0, 100.0, size=(10, 2)), index=years, columns=["Hopland", "Healdsburg"] ) # Two sites, each [time x realizations] df_hopland = pd.DataFrame(rng.gamma(2.0, 100.0, size=(10, 3)), index=years, columns=[0, 1, 2]) df_healdsburg = pd.DataFrame(rng.gamma(2.0, 100.0, size=(10, 3)), index=years, columns=[0, 1, 2])

In [4]:

from synhydro.core.ensemble import Ensemble, EnsembleMetadata

# Keyed by realization: each DataFrame is [time x sites]
data = {
    0: df_realization_0,
    1: df_realization_1,
}

# Or keyed by site: each DataFrame is [time x realizations]
data = {
    "Hopland": df_hopland,      # columns are realization labels
    "Healdsburg": df_healdsburg,
}

ensemble_custom = Ensemble(
    data,
    metadata=EnsembleMetadata(
        time_resolution="YS",
        description="Annual paleo reconstruction ensemble",
    ),
)
print(ensemble_custom)

from synhydro.core.ensemble import Ensemble, EnsembleMetadata # Keyed by realization: each DataFrame is [time x sites] data = { 0: df_realization_0, 1: df_realization_1, } # Or keyed by site: each DataFrame is [time x realizations] data = { "Hopland": df_hopland, # columns are realization labels "Healdsburg": df_healdsburg, } ensemble_custom = Ensemble( data, metadata=EnsembleMetadata( time_resolution="YS", description="Annual paleo reconstruction ensemble", ), ) print(ensemble_custom)

============================================================
Ensemble Summary
============================================================
Realizations: 3
Sites: 2
Time Period: 2000-01-01 to 2009-01-01
Description: Annual paleo reconstruction ensemble
============================================================

Note (Realization keys are integers): Realization identifiers are always integers. When you pass site-keyed data whose columns are integer-like strings (for example "1", "2" read from CSV headers), the labels are coerced to int automatically. Labels that cannot be coerced (for example "ens_a") are kept as-is and a warning is logged, because string realization keys are indistinguishable from site names if the dictionary is later used to build another Ensemble.

Set time_resolution in the metadata when constructing from raw data. Disaggregators check it to confirm the input frequency (for example "YS" for annual input to an annual-to-monthly NowakDisaggregator).

Summary statistics¶

summary reduces each site (or realization) to scalar statistics, and percentile returns time-varying quantiles across realizations:

In [5]:

stats = ensemble.summary(by="site")            # mean/std/min/max per site
bands = ensemble.percentile([10, 50, 90], by="site")
median_flow = bands["USGS-01434000"]["p50"]    # median across realizations
stats

stats = ensemble.summary(by="site") # mean/std/min/max per site bands = ensemble.percentile([10, 50, 90], by="site") median_flow = bands["USGS-01434000"]["p50"] # median across realizations stats

Out[5]:

	mean	std	min	max
site
USGS-01434000	4493.509654	2552.273730	356.231863	35034.024620
USGS-01438500	5119.105968	2891.242967	470.994406	30120.827090
USGS-01440000	101.255929	65.357896	4.748189	1694.498639
USGS-01463500	10551.470469	5731.316109	769.438852	68496.074252

Subsetting and resampling¶

Both return a new Ensemble and leave the original untouched:

In [6]:

subset = ensemble.subset(
    sites=["USGS-01434000", "USGS-01438500"],
    realizations=[0, 1, 2],
    start_date="2000-01-01",
    end_date="2010-12-31",
)

annual = ensemble.resample("YS")   # sums to the new frequency

print(subset)
print(annual)

subset = ensemble.subset( sites=["USGS-01434000", "USGS-01438500"], realizations=[0, 1, 2], start_date="2000-01-01", end_date="2010-12-31", ) annual = ensemble.resample("YS") # sums to the new frequency print(subset) print(annual)

============================================================
Ensemble Summary
============================================================
Realizations: 3
Sites: 2
Time Period: NaT to NaT
Generator: KirschGenerator
Description: Subset of ensemble
============================================================
============================================================
Ensemble Summary
============================================================
Realizations: 50
Sites: 4
Time Period: 2026-01-01 to 2055-01-01
Generator: KirschGenerator
Description: Resampled to YS
============================================================

Saving and loading¶

Ensembles serialize to HDF5. By default data is grouped by site (stored_by_node=True), which matches the layout expected by downstream tools such as Pywr-DRB:

In [7]:

ensemble.to_hdf5(OUT_DIR / "synthetic_flows.h5")

loaded = Ensemble.from_hdf5(OUT_DIR / "synthetic_flows.h5")
first_three = Ensemble.from_hdf5(
    OUT_DIR / "synthetic_flows.h5", realization_subset=[0, 1, 2]
)
print(loaded)
print(first_three)

ensemble.to_hdf5(OUT_DIR / "synthetic_flows.h5") loaded = Ensemble.from_hdf5(OUT_DIR / "synthetic_flows.h5") first_three = Ensemble.from_hdf5( OUT_DIR / "synthetic_flows.h5", realization_subset=[0, 1, 2] ) print(loaded) print(first_three)

============================================================
Ensemble Summary
============================================================
Realizations: 50
Sites: 4
Time Period: 2026-01-01 to 2055-12-01
Generator: KirschGenerator
============================================================
============================================================
Ensemble Summary
============================================================
Realizations: 3
Sites: 4
Time Period: 2026-01-01 to 2055-12-01
Generator: KirschGenerator
============================================================

Next steps¶

• Quantitative validation of an ensemble - Tutorial 05
• Plotting ensembles - Tutorial 06
• Full API reference - Core Data Structures

---

Previous: Disaggregator (Nowak) | Next: Verification & Validation