Layer 3 — Integration Examples¶
Cross-library patterns showing how ekalmX composes with the broader ecosystem.
With optax (EKP as optimizer)¶
EKI as a drop-in optax optimizer¶
import ekalmx
import optax
# EKI as drop-in optimizer — parallel to optax_bayes
optimizer = ekalmx.optax.eki(
forward_fn=simulator,
obs=observations,
noise_cov=R,
scheduler=ekalmx.DataMisfitController(),
n_ensemble=100,
)
# Compose with optax ecosystem
optimizer = optax.chain(
optimizer,
optax.clip_by_global_norm(1.0),
)
# Standard optax loop
state = optimizer.init(params)
for step in range(n_steps):
updates, state = optimizer.update(grads=None, state=state, params=params)
params = optax.apply_updates(params, updates)
With JAX/Equinox (Differentiable training)¶
Learning neural dynamics by backpropagating through the filter¶
import jax
import optax
import ekalmx
import equinox as eqx
# Neural ODE dynamics
class NeuralDynamics(ekalmx.AbstractDynamics):
net: eqx.nn.MLP
def __call__(self, state, t0, t1):
# Simple Euler step with learned RHS
dt = t1 - t0
return state + dt * self.net(state)
dynamics = NeuralDynamics(net=eqx.nn.MLP(in_size=40, out_size=40, ...))
# Loss: negative log-likelihood from filter
def loss_fn(dynamics, ensemble, obs_sequence, obs_op, R):
filter = ekalmx.filters.ETKF()
particles = ensemble
total_ll = 0.0
for obs_val, t0, t1 in obs_sequence:
# Forecast — differentiable through neural ODE
particles = eqx.filter_vmap(dynamics)(particles, t0, t1)
# Analysis — differentiable through Kalman update
result = filter.analysis(particles, obs_val, obs_op, R)
particles = result.particles
total_ll = total_ll + result.log_likelihood
return -total_ll
# Train with standard JAX/optax loop
optimizer = optax.adam(1e-3)
opt_state = optimizer.init(eqx.filter(dynamics, eqx.is_array))
@eqx.filter_jit
def train_step(dynamics, opt_state, ensemble, obs_sequence, obs_op, R):
loss, grads = eqx.filter_value_and_grad(loss_fn)(
dynamics, ensemble, obs_sequence, obs_op, R
)
updates, opt_state = optimizer.update(grads, opt_state)
dynamics = eqx.apply_updates(dynamics, updates)
return dynamics, opt_state, loss
With somax (Ocean state estimation)¶
Wrapping a somax shallow water model as ekalmX dynamics¶
import somax
import ekalmx
# somax forward model
grid = somax.ArakawaCGrid2D(nx=128, ny=128, dx=10e3, dy=10e3)
swm = somax.ShallowWaterModel(grid=grid, params=somax.SWMParams())
# Wrap as ekalmX dynamics
class OceanDynamics(ekalmx.AbstractDynamics):
model: somax.ShallowWaterModel
dt: float = eqx.field(static=True, default=300.0)
def __call__(self, state, t0, t1):
return self.model.integrate(state, t0=t0, t1=t1, dt=self.dt).ys[-1]
# Sparse observation operator (satellite tracks)
class SatelliteObs(ekalmx.AbstractObsOperator):
obs_indices: Int[Array, "N_y"]
def __call__(self, state):
# Extract SSH at satellite footprint locations
return state[self.obs_indices]
# Run LETKF
letkf = ekalmx.LETKF(
dynamics=OceanDynamics(model=swm),
obs_op=SatelliteObs(obs_indices=sat_indices),
localizer=ekalmx.GaspariCohn(radius=300e3),
inflator=ekalmx.RTPS(alpha=0.85),
config=ekalmx.FilterConfig(n_ensemble=40),
)
result = letkf.assimilate(init_ensemble, observations, t0=0.0, dt=3600.0)
With gaussx (Structured noise)¶
Diagonal and low-rank noise models via gaussx operators¶
import lineax as lx
import gaussx
import ekalmx
# Diagonal observation noise
R_diag = lx.DiagonalLinearOperator(jnp.full(n_obs, 0.01))
# Structured process noise (low-rank + diagonal)
Q = gaussx.operators.low_rank_plus_diag(
W=basis_vectors, # (N_x, r) — r dominant modes
d=jnp.full(n_x, 1e-4), # small isotropic baseline
)
# Both plug directly into ekalmX
result = filter.analysis(particles, obs, obs_op, obs_noise=R_diag)
With geo_toolz (Pre/post-processing pipeline)¶
Preprocessing satellite data and evaluating filter output¶
import geo_toolz as gt
import ekalmx
# --- Preprocessing ---
raw_obs = xr.open_dataset("satellite_ssh.nc")
pipeline = gt.Sequential([
gt.validation.harmonize_coords(),
gt.subset.bounding_box(lon=(-80, 0), lat=(20, 60)),
gt.regrid.linear(target_grid=model_grid),
gt.detrend.remove_climatology(climatology=ssh_clim),
])
obs_dataset = pipeline(raw_obs)
# Convert xarray → JAX arrays for ekalmX
obs_values = jnp.array(obs_dataset["ssh_anomaly"].values)
# --- Run ekalmX ---
result = letkf.assimilate(ensemble, obs_values, ...)
# --- Postprocessing ---
analysis_ds = xr.Dataset({"ssh": (["ens", "y", "x"], result.particles)})
metrics = gt.Sequential([
gt.metrics.rmse(reference=truth_dataset),
gt.spectral.psd(dim="x"),
])
evaluation = metrics(analysis_ds)
With xr_assimilate (xarray orchestration)¶
Using ekalmX as a compute backend for xarray-level DA¶
import xr_assimilate as xra
import ekalmx
# ekalmX provides array-level compute
filter_backend = ekalmx.LETKF(
dynamics=ocean_dynamics,
obs_op=satellite_obs,
localizer=ekalmx.GaspariCohn(radius=300e3),
config=ekalmx.FilterConfig(n_ensemble=40),
)
# xr_assimilate handles xarray bookkeeping
assimilator = xra.Assimilator(
backend=filter_backend,
observations=obs_dataset,
background=background_dataset,
output_dir="./analysis/",
)
result_ds = assimilator.run() # xr.Dataset with analysis fields
Tutorial Index¶
Full end-to-end tutorials (in notebooks/):
| Tutorial | Layer | Topic |
|---|---|---|
| 01_primitives | L0 | Ensemble statistics, Kalman gain, localization from scratch |
| 02_linear_gaussian | L1 | EnKF on a linear Gaussian system — verify against analytic KF |
| 03_lorenz96_etkf | L1-L2 | ETKF on Lorenz-96 with localization and inflation |
| 04_eki_calibration | L1-L2 | Calibrate L96 forcing parameter F via EKI |
| 05_differentiable_enkf | L1 | Learn neural dynamics by backpropagating through EnKF |
| 06_ocean_letkf | L2 | LETKF with somax shallow water model |
| 07_eki_optax | L2 | EKI as optax transform, composed with schedules |
| 08_smoother | L1-L2 | EnKS backward pass, compare filter vs smoother |