prophys 0.4.0

Risk models you can differentiate.

Overview

Write the geometry, the physics and the uncertainty of a problem as one Python model. prophys compiles it into a single JAX graph. You sample it, fit it to measurements and optimise a design with the same gradient, and every number comes back with its error bar.

A pipeline route, an abatement fraction, a copula correlation and the weights of a learned surrogate are all parameters of that one graph, so a single jax.grad pass moves all of them. NumPyro and Distrax give you distributions. prophys builds the models around them: geometry, fields, physics and learned components wired together, with units checked where they meet.

It combines:

• A symbolic expression core (Expr, Param, Input, RandomVariable) with NumPy-style operator overloading. Models are built with ordinary arithmetic and compile to pure jax.numpy.

• Differentiable geometry structures in shared coordinate frames: PointList, LineList, Polygon/PolygonList (signed distance, containment, area), Polyhedron (3D halfspace volumes), raster Fields with bilinear interpolation, and evaluation Grids.

• Transformations: linear/affine/polynomial, smoothable PiecewiseLinear, Decay, Logistic dose-response, and interpolated TableLookup1D/2D characteristic curves.

• Distributions: 12+ univariate families, circular (VonMises, WrappedGaussian), multivariate Gaussian (Cholesky-parameterized), mixtures, empirical/KDE, and Gaussian copulas with always-valid learnable correlation.

• Probabilistic inputs anywhere in the graph: a RandomVariable (e.g. wind speed/direction) can feed the physics, not just the output noise. The compiled model marginalizes it by reparameterized Monte Carlo with gradients intact.

• A training engine (Optax): calibrate one attribute, finetune the whole model jointly across attributes with shared parameters, fit_distribution for standalone/multivariate distribution learning, with trainable/frozen parameter control.

• Design optimization with penalty constraints, where bounded parameters are constrained structurally (sigmoid/softplus bijectors) rather than by projection.

• An observation model: left- and right-censored, interval-censored and missing records each contribute their correct likelihood term instead of a substituted number, and prophys.io reads them out of CSV, wind roses, GeoTIFF and NetCDF with units converted at the boundary.

• Uncertainty on every number: Monte-Carlo standard errors on each sampled statistic (batch means, Maritz-Jarrett quantiles, influence-function CVaR), and observed-information standard errors, confidence intervals, profile likelihoods and delta-method propagation on each fitted parameter.

• Sensitivity analysis: exact local derivatives and dimensionless elasticities by autodiff, and global Sobol variance decomposition separating direct effects from interactions.

• Model validation: AIC/AICc/BIC, WAIC and PSIS-LOO with Pareto-k diagnostics, model comparison with the standard error of the difference, PIT calibration, CRPS, and posterior predictive checks.

• Posterior inference through NumPyro: NUTS over the compiled model, with split-r-hat, effective sample sizes and divergence reporting.

• Units that are checked: dimensions, SI prefixes, compound expressions and affine temperature scales, enforced at frame and attribute boundaries at model-build time.

• A standardized export format (ModelPackage): quantized PMF, histogram, or JAX-free sampler representations, JSON round-trippable, consumer-agnostic, carrying a run manifest that records what produced the result.

• Learned surrogates (prophys.nn): networks, differentiable tree ensembles, field transformers and the DPOT neural operator, pretrained on their own data and fine tuned inside the model with a feedback dial between cut and joint estimation. Monotone, convex, Lipschitz and bounded layers, physics layers written as prophys expressions that compile into the network and export with it, physics informed residuals with exact Taylor mode derivatives, quantile, Gaussian and mixture heads, Laplace and anchored ensemble uncertainty that enters the Monte Carlo draw, extrapolation guards, XGBoost, LightGBM, Transolver and DPOT importers, and a prior fitted network that predicts from a few dozen rows without fitting. All numerics run in optimized native kernels at TensorFlow speed, DPOT trains faster than its PyTorch reference, and results are identical on CPU and GPU.

• Explainability (prophys.explain): exact TreeSHAP, exact and kernel Shapley values with interactions, expected and integrated gradients, and Shapley attribution of variance and tail risk through a whole model.

• Diagnostics with stable codes, locations and hints for layer dimensions, training convergence, licence limits, symbolic surrogates and physics primitives. The convergence check also watches the parameters, so a run that stopped on a plateau before a later drop is reported as unfinished.

• A runtime library (prophys-runtime, a shared library with a C interface) that evaluates exported surrogates without Python, bit identical to prophys.

Installation

pip install prophys

Installs ready to run on Linux, macOS and Windows (x86_64 and arm64) with its optimized native kernels included; there is nothing to compile.

Free tier / licensing: models declaring up to 250 structural objects in total (points, polygon vertices, segments, raster cells, ...) run without any license. This does not count Monte-Carlo sample counts, observation counts, or training iterations, only declared model structure. Larger models require a signed license (set via the PROPHYS_LICENSE environment variable). Contact RhineQC GmbH to obtain one.

With optional extras:

pip install "prophys[viz]"        # matplotlib plotting helpers (prophys.plot)
pip install "prophys[geo]"        # shapely/pyproj geodata import
pip install "prophys[inference]"  # NumPyro bridge

Requirements: Python ≥ 3.11 (JAX and Optax install automatically as dependencies). Python 3.10 is not supported: jaxlib has not published a cp310 wheel since 0.9.0.

Quick Start

import jax
import jax.numpy as jnp
import optax
import prophys as prp

# Geometry in a shared frame: the pipeline route (optimizable) and a receptor
site = prp.Frame("site", units="m")
route = prp.Polyline(
    prp.Param("route", shape=(4, 2), init=jnp.array([[0.0, 0.0], [300.0, 60.0], [650.0, -30.0], [900.0, 0.0]])),
    site,
)
house = prp.Point(jnp.array([500.0, 250.0]), site)

# Physics: a leak can occur anywhere along the route, so its location is a
# continuous, reparameterized random variable pushed through the route's
# differentiable point_at(t), not a hand-picked set of candidate points
leak_frac = prp.RandomVariable("leak_frac", prp.Uniform(low=0.0, high=1.0))
leak_point = route.point_at(leak_frac)

# Gas dispersion from the leak to the house, then concentration -> a
# probability of health issues (dose-response curve)
distance = prp.norm(house.coords[0] - leak_point)
concentration = prp.Param("emission_rate", init=1.5) * prp.exp(-distance / 150.0)
health_risk = prp.PiecewiseLinear([0.0, 0.05, 0.2, 0.5], [0.0, 0.05, 0.3, 0.85], smooth=0.02)(concentration)

# Probabilistic attribute + standardized model object
incident = prp.UncertainAttribute("health_incident", prp.Bernoulli(prob=prp.clip(health_risk, 1e-4, 1 - 1e-4)))
model = prp.ProbabilityModel(incident)
compiled = model.compile(mode="opt", n_samples=200)

# Probability of a health incident, marginalized over every leak location
# along the route, with a gradient with respect to every vertex of the route
p_incident = compiled.expectation("health_incident")
grad_wrt_route = jax.grad(
    lambda params: compiled.expectation("health_incident", params=params)
)(compiled.default_params())

# Calibrate against observations ...
result = prp.calibrate(compiled, "health_incident", observations, n_steps=300)

# ... or optimize the route itself, directly against the marginalized
# probability, with plain JAX/optax (prp.optimize is only a convenience
# wrapper around the same pattern, not required)
opt = optax.adam(0.05)
raw_params = compiled.default_params()
opt_state = opt.init(raw_params)

for _ in range(200):
    loss, grads = jax.value_and_grad(
        lambda params: compiled.expectation("health_incident", params=params)
    )(raw_params)
    updates, opt_state = opt.update(grads, opt_state)
    raw_params = optax.apply_updates(raw_params, updates)

Uncertainty at the inputs, not just the output

wind_speed = prp.RandomVariable("v", prp.Weibull(scale=6.0, concentration=2.0))
wind_dir   = prp.RandomVariable("phi", prp.VonMises(loc=prp.deg2rad(35.0), kappa=3.0))

# Any expression downstream of these is a random quantity. The compiled
# model marginalizes them via reparameterized Monte Carlo, inside one
# JAX trace, so gradients flow back through the marginalization.

Training the model

# Fit one attribute
prp.calibrate(compiled, "exposure", observations)

# Jointly train the entire model. Shared parameters pool information
# across every attribute's observations
prp.finetune(compiled, {"exposure": obs_e, "complaints": obs_c}, frozen=["emission_rate"])

# Learn a standalone (also multivariate) distribution from data
mvn = prp.MultivariateGaussian(
    mean=prp.Param("mu", shape=(2,), init=jnp.zeros(2)),
    cholesky=prp.Param("L", shape=(2, 2), init=jnp.eye(2)),
)
prp.fit_distribution(mvn, data)

Records that are not exact numbers

A concentration below a detection limit, a specimen still unfailed at the end of a test, a reading binned at recording, a failed sensor: each carries a different likelihood term, and substituting a number for any of them changes the fitted parameters.

obs = prp.Observations.detection_limit(concentrations, limit=0.5)   # "< 0.5"
obs = prp.Observations.right_censored(times, censored=still_running)
obs = prp.Observations.interval(lower_edges, lower_edges + width)
obs = prp.Observations.from_arrays(values, missing=np.isnan(values))

prp.finetune(compiled, {"c": obs})    # accepted anywhere an array was

prophys.io reads them straight out of the formats they arrive in, such as a CSV with <0.5 and ND entries, a wind rose in degrees, a GeoTIFF or a NetCDF extract, with units declared per column and converted at the boundary.

The uncertainty of the numbers

Every reported statistic comes with its own error, and every fitted parameter with its own interval.

compiled.assess("loss", "cvar", level=0.99)   # estimate + Monte-Carlo error
prp.required_samples(estimate, 0.01)          # draws needed for 1% precision

u = prp.parameter_uncertainty(compiled, {"y": data}, result=result)
u.interval("sigma")                           # standard errors + intervals
u.collinear_pairs()                           # what the data cannot separate
prp.propagate(u, lambda p: p["mu"] + 1.645 * p["sigma"])   # into any quantity

The Monte-Carlo error uses batch means for a mean, the Maritz-Jarrett order-statistic estimator for a quantile, and the influence function for a CVaR, the last one including the uncertainty of the VaR threshold itself. Parameter errors come from the observed information at the optimum, mapped into constrained space so a positive scale never gets an interval straddling zero; an unidentified direction is reported as such rather than as numbers.

What the answer depends on

prp.local_sensitivity(compiled, "power")      # exact derivatives, elasticities
prp.sobol_sensitivity(compiled, "power")      # variance shares + interactions
compiled.conditional_expectation("power", {"wind": np.array([4., 8., 12.])})

Whether the model is any good

prp.information_criteria(compiled, {"y": data}, result=result)   # AIC/AICc/BIC
prp.waic(compiled, "y", data, draws)                             # predictive
prp.loo(compiled, "y", data, draws)                              # + Pareto k
prp.compare({"simple": a, "richer": b})                          # with SE of Δ
prp.calibration(compiled, "y", data)      # under-dispersed? biased?
prp.crps(compiled, "y", data)             # the whole distribution, in its units
prp.posterior_predictive_check(compiled, "y", data, statistic="skew")

prp.report_card(compiled, "y", data, result=result)   # all of the above

The posterior, not just the mode

from prophys.inference import sample_posterior     # pip install "prophys[inference]"

posterior = sample_posterior(compiled, {"y": data}, init_params=result.raw_params)
print(posterior.summary())        # credible intervals, r_hat, ESS, divergences
prp.waic(compiled, "y", data, posterior.thin(400))

Units that are checked

site = prp.Frame("site", units="m")
plan = prp.Frame("plan", units="km")
prp.FrameTransform(plan, site).apply(pts)   # converts; it is not the identity

power = prp.UncertainAttribute("power", dist, unit="kW")
power.convert(1000.0, "MW")   # 1.0
power.convert(1.0, "kWh")     # UnitError: power is not energy

Export

pkg = prp.export.export_model(
    compiled,
    params=result.raw_params,
    observations={"y": data},
    diagnostics=prp.report_card(compiled, "y", data, result=result),
)
pkg.save("model.json")       # versioned, consumer-agnostic interchange format

Alongside the quantized PMF / histogram / JAX-free sampler, format 1.1 carries a run manifest: the prophys version, a structure hash of the model, the parameters in both constrained and unconstrained form, the seed and Monte-Carlo settings, the licensing regime, the observation counts by kind, and any attached diagnostics. A later run can verify against it:

loaded = prp.export.ModelPackage.load("model.json")
loaded.verify(model.compile())                     # raises if the model changed
compiled = prp.export.reproduce(loaded.manifest, model)

Downstream consumers implement their own import against the documented ModelPackage shape. prophys carries no consumer-specific dependencies.

Examples

Runnable end-to-end models in examples/, each documented step by step and running on the free tier:

• gas_dispersion/: pipeline-leak health risk, a Gaussian plume from a leak anywhere along the route, calibration and abatement by gradient.
• runway_crosswind/: crosswind and wind-shift operability of two runways under a coupled, circular wind process.
• electric_grid/: fault cascades in a radial feeder, line reinforcement and substation siting optimised together.
• power_flow/: transmission-grid security, PTDF and LODF from the DC power flow, N-1 screening, cascading overloads and line uprates.
• wind_farm/: blade element momentum rotor, directional wind climate, wake losses, AEP with P50/P90 and layout optimisation.
• airfoil_loads/: gust loads and stall risk of a light aircraft with a calibrated lift polar.
• rail_fatigue/: rainflow fatigue of a railway bridge detail, calibrated to strain gauges, with reliability over the service life.
• groundwater_wells/: Theis drawdown of a well field beside a protected fen, calibrated to a pumping test, with pumping rates under a chance constraint.
• pharmacokinetics/: population pharmacokinetics from sparse, censored samples and dose-regimen design.
• insurance_layers/: exact aggregate losses by Panjer recursion, stop-loss layer pricing and an optimal retention.
• learned_surrogate/: a physics-informed network surrogate of a battery thermal model inside a pack risk model.

Documentation

Full docs are at docs.rhineqc.com/distribution/prophys, covering concepts (symbolic graphs, frames, eval/opt modes, marginalization), the full structure/distribution/transformation catalogs with plots, units and dimensional checking, reading measurement data, calibration and design optimization, Monte-Carlo and parameter uncertainty, sensitivity analysis, model validation, posterior inference, step-by-step worked examples, and the export format specification.

License

prophys is proprietary software by RhineQC GmbH. It is intended for internal or explicitly licensed distribution.