PyCausalSim 0.3.0: Causal Simulation for Monthly Forecasts, Resampled History and Honest Backtests

An open-source Python library for causal inference and simulation, now with whole-episode resampling, window claims and rolling backtests, built and tested on a real natural gas model.

Note: I am the creator of PyCausalSim. You can find it on GitHub at github.com/Bodhi8/pycausalsim.

PyCausalSim 0.3.0 is out. It is free, open source under the MIT license, and one command away, installed straight from GitHub:

pip install git+https://github.com/Bodhi8/pycausalsim.git

PyCausalSim answers two kinds of question with the same machinery. Looking backward, it asks what caused an outcome already in your data. Looking forward, its futures module asks what happens under a scenario nobody has seen yet: you state the causal graph, unroll it over a horizon, run Monte Carlo draws over the uncertain parts, and get ranges, scenario comparisons and dated claims you can score later.

This release came out of a real build. For my article on improving LNG and natural gas forecasting, I modeled a Henry Hub winter month by month, with the NYMEX futures curve as the center and storage surprises as the causal channel. Every rough edge I hit became a feature or a fix. The result is a small release with a clear purpose: make forward-looking causal models easier to build correctly, and harder to get subtly wrong.

What’s in 0.3.0

New in pycausalsim.futures:

• empirical_path() resamples whole historical episodes, such as past winters, so the pattern across months stays intact. The pick is shared across scenarios on the same draw, and optional weights let outside probabilities steer it.
• parameter(..., "choice") gives you finite priors: a set of values with optional probabilities.
• Scenario.path() sets a lever month by month in one call.
• RunResult.window() returns each draw's average, total, minimum or maximum over a span, and claims can test it directly.
• rolling_backtest() scores a model from many forecast origins and pools the result, with each origin's own scorecard kept.
• Monthly horizons are first-class: steps can be any sortable label, such as 202611 for November 2026, and c.step and scm.steps name them that way.
• A new worked example, examples/henry_hub_winter.py, ships with EIA data and a regression test.

Fixes and housekeeping:

• Uplift modeling’s AUUC and Qini metrics no longer crash on NumPy 2, which removed the function they used.
• Python 3.10 or later, matching what CI actually tests, with one version number in one place.
• The README now lists only the causal discovery methods that exist.
• CI runs every example on every push, so a broken example cannot ship again.
• The test suite grew from 59 to 72 tests.

Why these features: what a real monthly model needed

Forecasting a natural gas winter exposed four places where the library made the right thing harder than it should be.

Weather comes in whole winters. A cold December usually means an expensive January and February too. If each month’s weather is drawn independently, the winter-average range comes out far too narrow and every risk number built on it is wrong. Version 0.2 had no direct way to replay a past winter as a unit, so I faked it with a uniform draw rounded to an index. empirical_path() does it properly: each draw replays one historical episode end to end, and every scenario in a comparison sees the same episode in the same draw, so differences between scenarios are paired.

The weights are the point. empirical_path() takes optional weights over the past episodes. That is how outside probabilities enter a causal simulation without retraining anything: give more weight to the past winters that look like this year's seasonal outlook, and the price ranges, scenarios and claims shift with it. In the gas work, a classifier built from prices added nothing, but a forecast from outside the market, such as a weather outlook, is exactly the kind of input that can.

Decisions are about spans, not single months. A fuel budget cares about the November-to-March average. In 0.2, a claim could only test one step, so I carried a running sum as an extra state variable. RunResult.window() returns each draw's aggregate over any span, and claim tests can call it: "the winter average settles below the strip" is now one line.

Honest backtests need many origins. The test that matters for a seasonal model is retrodiction from each past October, with that winter left out of the calibration. I wrote that loop by hand. rolling_backtest() now pools any set of origin-specific models and reports pooled error and coverage, plus each origin's own scorecard.

One lesson did not need code, only documentation. While backtesting hedging decisions I pulled each month from a separate simulation run. Each run draws fresh parameters, so the months stopped moving together and the ranges collapsed. The fix was one line, but the mistake was easy to make, so the README now says it plainly: take every quantity for one analysis from one run. compare() and window() keep you there.

Quickstart

A toy winter model that uses every new piece and needs no data files:

import numpy as np
from pycausalsim.futures import TemporalSCM, Scenario, compare

months = [202611, 202612, 202701, 202702, 202703]

# Five past winters of storage surprises (negative = tighter), one row per winter, one column per month.
past_winters = np.array([
    [ 10,  -40,  -60, -20,  0],
    [ 20,   30,   10,  40, 30],
    [-30,  -90, -120, -50, -10],
    [  0,   10,  -20,  60, 20],
    [ 15,  -20,   40, -30, 10],
])

scm = TemporalSCM(horizon=months, step="month", draws=10000, seed=7)
scm.parameter("slope", "normal", loc=0.002, scale=0.0005)        # price response per unit
scm.parameter("regime", "choice", values=[0.9, 1.0, 1.2], p=[0.3, 0.5, 0.2])
scm.lever("extra_demand", default=0.0)

# Replay whole winters; give the tight 3rd winter double weight, as a cold outlook might.
scm.empirical_path("shock", past_winters, weights=[1, 1, 2, 1, 1])

scm.state("tightness", lambda c: c.lag("tightness", 1) - c["shock"] + c["extra_demand"],
          parents=["tightness:lag(1)", "shock", "extra_demand"], init=lambda p: 0.0)
scm.add_equation("price", lambda c: 3.0 * c.p.regime * np.exp(c.p.slope * c["tightness"]),
                 parents=["tightness"])

scm.claim("winter_above_3", "Winter average price above $3.00",
          test=lambda r: r.window("price", 202611, 202703) > 3.0, resolve_by=202703)

base = Scenario("Base")
surge = Scenario("Demand surge").path("extra_demand", {202612: 20, 202702: 0})
cmp = compare(scm, [base, surge], target="price")

for name in cmp.runs:
    w = cmp[name].window("price", 202611, 202703)
    print(f"{name:13s} winter average {w.mean():.2f}  90% [{np.percentile(w, 5):.2f}, {np.percentile(w, 95):.2f}]")
print(cmp.falsifiable_claims(scenario="Base"))

It prints:

Base          winter average 3.53  90% [2.47, 5.34]
Demand surge  winter average 3.75  90% [2.60, 5.76]
[202703] p=0.59  Winter average price above $3.00  (Base)

The demand surge runs only in December and January, through Scenario.path, but its effect carries into February and March through the state variable, which is exactly the kind of lagged consequence a single-period model misses.

The worked example: Henry Hub, winter 2026–27

examples/henry_hub_winter.py is the full model behind my LNG forecasting article, with public EIA data bundled so it runs out of the box:

python examples/henry_hub_winter.py

It anchors on the NYMEX futures strip from October 23, 2026, estimates how much a storage surprise moves price relative to that curve, resamples weather from the last twelve winters with empirical_path, and runs scenarios on the same draws. What it reports:

• Market view: a November-to-March average of $2.70, with a 90% range of $1.82 to $3.82, against a strip of $2.97. A 70% probability that the winter settles below the strip.
• January: a mean of $3.22 and a 90% range of $1.97 to $4.92.
• LNG feedgas 1.5 Bcf/d above expectations all winter: a winter average of $3.01.
• A replay of the 2013–14 polar vortex: $5.06, flagged by PyCausalSim as extrapolated, because that winter’s storage path lies outside the range the price response was estimated on.
• A hedging decision: to keep the 95th-percentile winter cost within 10% of the strip, lock 65% of volume, for an expected cost of $2.88.

The example also runs winter_backtest(), a leave-one-winter-out test from each October curve from 2014-15 to 2023-24, through rolling_backtest. Across 40 winter months, the average miss was 22.3% and 88% of realized prices fell inside the 90% range, close to where a calibrated model should be. The test suite pins these numbers, so a change that moves them fails CI.

One honest note from the article carries over: for a single-number forecast, most of the model’s edge over the curve came from correcting the curve’s winter premium, and a simple premium shift did about as well. The simulation earns its keep elsewhere, in calibrated ranges, priced scenarios and decisions tied to an explicit risk limit. The full walkthrough is in How to Use PyCausalSim to Improve LNG and Natural Gas Forecasting.

What’s next

The next step for the gas model is also the next step for the library: feeding NOAA’s seasonal temperature outlook into empirical_path weights, so a published probability forecast steers the simulation directly. After that, the roadmap is:

• An agent layer for futures, so heterogeneous households, firms or regions can sit inside a temporal causal model instead of beside it.
• Geography, starting with regional prices and flows between regions.
• Uncertainty over structure, not only coefficients, by sampling over a set of plausible graphs.
• Scenarios in YAML, so people who do not write Python can author and review them.
• FCI for causal discovery with latent confounders, the one method the README used to claim and the code did not have.

Try it

pip install git+https://github.com/Bodhi8/pycausalsim.git
python -c "import pycausalsim; print(pycausalsim.__version__)"

The code, examples and changelog are at github.com/Bodhi8/pycausalsim. If you build a forward-looking model with it, in energy, healthcare, marketing or anything else, open an issue or a discussion and tell me what broke. This release exists because a real model broke a few things, and the next one will too.

About the Author

Brian James Curry is a Kansas City AI researcher, applied scientist and data scientist with more than twenty years in marketing modeling and measurement. His work today spans marketing mix modeling, causal inference, machine learning and healthcare intelligence. He is the founder of Vector1 Labs, an independent applied AI research lab. His marketing measurement work covers marketing mix modeling, incrementality testing and budget allocation, and his Causal Care series applies causal machine learning to health system decisions: which programs work, for whom, and what they are worth. He also founded KC AI Lab, which began in 2016 as the Kansas City Machine Learning Group and has grown to more than 1,500 members.

His career includes leading analytics, data science and machine learning at Hallmark, Garmin, Koch Industries, Tractor Supply Company and Vail Resorts. His open-source research tools on GitHub include Papilon (marketing mix modeling and complex systems), PyCausalSim (causal discovery through simulation) and tincast (forecasting for AI-critical commodities), along with MeaningFlow and Daedalus. A Kansas State University alum from a small Kansas farming town, he writes about marketing mix modeling, causal inference, agentic AI, healthcare AI and the economics of AI on Medium and Substack. Contact: brian@vector1.ai.