One question, carried all the way · step 1 of 8
The question, written down
What are we actually deciding, and what would 'yes' have to beat?
A logistics operator runs 40 depots in one region and 500 nationally. Should the standing weekly maintenance schedule go from nothing to 60 hours per depot? Before any data is touched, the decision has to be written down as a quantity — because the number that answers it is not the number an experiment can most easily measure.
Two estimands, not one. The experiment can measure a first-week lift on individual depots; the decision turns on a steady-state weekly lift across the region. They are different windows and different levels, and keeping them apart from the first line is what stops the report answering the easy question and calling it the hard one.
surface_world is a simulator, and it hands back the spec for free so this story can get moving. With your own data there is no simulator and the spec is the first thing you write -- the porting tutorial does exactly that, from a dataframe, with nothing hidden.
from axiom.core import Intervention, Population, TimeWindow
from axiom.estimands import Estimand, Level, Quantity
from axiom.sim import DosePlan, surface_world
from axiom.surface import GeometricCarryover, HillKernel
DOSE, NO_DOSE = 60.0, 0.0 # maintenance hours per depot-week
VALUE_PER_POINT, COST_PER_HOUR = 400.0, 45.0
REGION_DEPOTS, MAX_LAG, SEED = 40, 4, 0
world = surface_world(
n_units=REGION_DEPOTS, n_periods=52, treatments=("a",),
kernels=HillKernel(reference_dose=50.0, amplitude_scale=10.0),
carryover={"a": GeometricCarryover(max_lag=MAX_LAG)},
doses=DosePlan(scale=50.0, spread=0.8, zero_fraction=0.05),
intercept="shared",
truth={"beta_a": 10.0, "alpha": 5.0, "k_a": 50.0, "s_a": 2.0, "lam_a": 0.5},
noise_sd=2.0, seed=SEED,
)
spec = world.spec
def contrast(name, window, level):
return Estimand(
name=name, quantity=Quantity(kind="contrast"),
treatment=spec.treatment("a"),
intervention=Intervention(doses={"a": DOSE}),
reference=Intervention(doses={"a": NO_DOSE}),
outcome=spec.outcome, population=Population(name="region_depots"),
window=window, level=level, dimension=spec.outcome_dimension,
)
experiment_estimand = contrast(
"first_week_lift",
TimeWindow(start=0, stop=1, basis="cumulative"),
Level(unit="individual"),
)
decision_estimand = contrast(
"steady_state_weekly_lift",
TimeWindow(start=MAX_LAG, stop=52, basis="per_period"),
Level(unit="cluster"),
)
BREAK_EVEN = DOSE * COST_PER_HOUR / VALUE_PER_POINT
print(f"the experiment can measure : {experiment_estimand.name}")
print(f"the decision turns on : {decision_estimand.name}")
print()
print(f"{DOSE:.0f} hours at {COST_PER_HOUR:.0f} USD "
f"= {DOSE * COST_PER_HOUR:,.0f} USD per depot-week")
print(f"break-even steady-state lift: {BREAK_EVEN:.2f} index points per depot-week")
the experiment can measure : first_week_lift the decision turns on : steady_state_weekly_lift 60 hours at 45 USD = 2,700 USD per depot-week break-even steady-state lift: 6.75 index points per depot-week
Every line above was captured by running this code at build time. The steps share one namespace, so the code on this page is the code you would type after the steps before it.