Ness, Causal AI, in StatsPAI¶
Robert Osazuwa Ness's book (Manning, 2025) teaches the graph side of causal inference: draw the DAG, test it against data, fit a model on it, intervene on the model, reason about counterfactuals, identify an estimand, and only then estimate. It works in pgmpy, y0, DoWhy and Pyro. This guide maps the chapters to the StatsPAI calls that give the same numbers and says where StatsPAI makes a different choice.
The notebooks and data are at https://github.com/altdeep/causalML. They
are not redistributed with StatsPAI. The examples below read the CSV files
from its datasets/ folder.
Chapter by chapter¶
| chapter | the book computes | StatsPAI |
|---|---|---|
| 3 | a DAG as a list of edges | sp.dag("A -> E; S -> E; ..."); a Graphviz digraph { ... } string is read as is |
| 3 | causal Markov kernels from data (pgmpy fit) |
sp.bayes_net(dag, df); .cpt(node) |
| 3 | the same with a Dirichlet prior | sp.bayes_net(dag, df, prior=1) |
| 3 | a conditional query by variable elimination | net.query("E", evidence={"T": "train"}) |
| 4 | d-separation | g.d_separated(x, y, {...}) |
| 4 | every independence the graph implies | g.implied_independencies() |
| 4 | chi-square tests of those on categorical data | g.test_implications(df) (test="g-test" for the likelihood ratio) |
| 6 | a structural model with assignment functions | sp.bayes_net(dag, cpts={node: function}) for discrete models; sp.SCM for continuous ones |
| 7 | the effect an experiment would show, from observational data | g.adjustment_sets(x, y), then sp.aipw / sp.ipw / sp.g_computation(by_arm=True) |
| 7 | graph mutilation, do |
g.do("X"); net.do(X="high"); net.query(..., do={...}) |
| 9 | counterfactuals on a parallel-world graph | net.counterfactual(variables, evidence=, do=) |
| 9 | abduction, action, prediction with continuous noise | sp.SCM(...).counterfactual(evidence, intervention) |
| 10 | is P(Y \| do(X)) identified, and by what formula |
sp.identify(g, x, y): .identifiable, .estimand |
| 10 | the estimand evaluated on data | sp.identify(g, x, y).estimate(df) |
| 10 | a counterfactual query: the effect of treatment on the treated | sp.identify_counterfactual(g, event, given=) |
| 11 | backdoor, front-door and instrumental strategies | g.recommend_estimator(x, y) lists all the graph licenses |
| 11 | regression, stratification, matching, weighting | sp.regress, sp.match, sp.ipw, sp.aipw |
| 11 | double machine learning, T-learner | sp.dml, sp.metalearner(learner="t") |
| 11 | front door by two regressions | sp.front_door(..., outcome_model="additive") |
| 11 | instrument, regression discontinuity | sp.iv, sp.rdrobust |
| 11 | refutation tests | sp.refute(estimator, df, y=, treat=, covariates=, method=) |
| 11 | sensitivity to an unobserved confounder | sp.sensemakr, sp.evalue |
| 12 | expected utility of an action, seen and done | net.expectation("U", evidence={...}) and net.expectation("U", do={...}) |
Chapters 2, 5 and 13, the flow-based models of chapter 6, the image examples of chapter 9 and the Bayesian neural model of chapter 11 are deep generative models in Pyro. StatsPAI does not cover them.
Build the graph, then test it¶
import pandas as pd
import statspai as sp
df = pd.read_csv("datasets/transportation_survey.csv") # six categorical columns
g = sp.dag("A -> E; S -> E; E -> O; E -> R; O -> T; R -> T")
g.implied_independencies() # eleven statements, the list dagitty gives
tests = g.test_implications(df) # chi-square for labels, Fisher z for numbers
tests.sort_values("p_value").head()
Each row is one independence the graph commits to. statistic, df and
p_value equal dagitty::localTests(type = "cis.chisq"). p_holm adjusts
for testing them all. cramers_v is there because, as the book shows in
listings 4.8 and 4.9, a p-value shrinks with the sample while the size of
the departure does not: with enough rows every such test rejects a graph
that is only approximately right.
Fit the graph and intervene on it¶
net = sp.bayes_net(g, df)
net.cpt("T") # P(T | O, R)
net.query("E", evidence={"T": "train"}) # seeing
net.query("T", do={"E": "uni"}) # doing
A table can also be written by hand, or given as a function of the parents. That is how a decision problem or a structural model is entered:
net = sp.bayes_net(
"C -> X; C -> Y; X -> Y; Y -> U",
cpts={
"C": {"bear": 0.5, "bull": 0.5},
"X": {"bear": {"debt": 0.8, "equity": 0.2},
"bull": {"debt": 0.2, "equity": 0.8}},
"Y": {("bear", "debt"): {"failure": 0.3, "success": 0.7},
("bull", "debt"): {"failure": 0.9, "success": 0.1},
("bear", "equity"): {"failure": 0.7, "success": 0.3},
("bull", "equity"): {"failure": 0.6, "success": 0.4}},
"U": lambda Y: -1000 if Y == "failure" else 99000,
},
)
net.expectation("U", evidence={"X": "debt"}) # 57000: what debt-financed firms earn
net.expectation("U", do={"X": "debt"}) # 39000: what choosing debt earns
Keys with several parents are tuples in the alphabetical order of the
parents' names (("C", "X") above).
Two things the network tells you that a table of numbers would not:
- A parent configuration that never occurs in the data has no frequencies.
A query whose answer depends on one warns that positivity fails, and
net.unseencounts them.prior=1makes the smoothing explicit. net.counterfactual(...)is refused unless every node other than the roots is a deterministic function of its parents. Tables alone do not say how one unit's value would change. Move each node's noise into a root and give the node a function, as the book does for Monty Hall.
Identify, then evaluate¶
When a confounder is unobserved the network cannot be fitted. The estimand can still be evaluated if the effect is identified:
g = sp.dag("X -> M -> Y; U -> X; U -> Y", latent=["U"])
res = sp.identify(g, "X", "Y")
res.estimand # 'sum_{M} [P(M | X) * sum_{X'} [P(X') * P(Y | M, X')]]'
res.estimate(df) # P(Y | do(X)) from the observed columns only
Rows come back missing, with a warning, for treatment values whose formula needs a conditional probability that has no data behind it.
Counterfactual queries¶
The effect of treatment on the treated compares what the treated got with
what they would have got untreated, P(Y_{X=0} = y | X = 1). The evidence
and the outcome live in different worlds, so sp.identify does not apply.
g = sp.dag("T -> W -> A; B -> V -> A; C -> T; C -> A; C -> B") # listing 10.7
res = sp.identify_counterfactual(g, [("A", 1, {"T": 0})], given=[("T", 1)])
print(res.summary())
res.estimate(df) # when res.from_observational_data
An event is (variable, value) for the world as it is and (variable,
value, {intervened: value}) for a world in which something was set. The
result gives the formula in interventional distributions, says for each
term whether observational data can replace the experiment, and evaluates
it.
Some counterfactuals no experiment can settle. The probability of
necessity, P(Y_{X=0} = 0 | X = 1, Y = 1), needs the joint behaviour of
Y in two worlds:
sp.identify_counterfactual(sp.dag("X -> Y"), [("Y", 0, {"X": 0})],
given=[("X", 1), ("Y", 1)]).identifiable # False
For those, state the structural model (sp.bayes_net(...).counterfactual,
sp.SCM) or bound the quantity.
One graph, three strategies¶
The book's online-game example (chapter 11) is identified three ways at once. Ask the graph:
g = sp.dag(causal_graph, latent=["Prior Experience"]) # the book's DOT string
rec = g.recommend_estimator("Side-quest Engagement", "In-game Purchases")
print(rec.summary())
Column names with spaces or hyphens go in backticks or Q("..."):
sp.regress("`In-game Purchases` ~ `Side-quest Engagement`", data=df).
The recommendation is the backdoor adjustment, and the alternatives list
the front door through Won Items and the instrument, each with a call
that runs. The strategies rest on different assumptions, so agreement among
them is evidence for the graph and disagreement is evidence against it. On
the book's data:
| strategy | call | estimate |
|---|---|---|
| backdoor, regression | the recommended sp.regress(...) |
178.09 |
| backdoor, doubly robust | sp.aipw(df, y=, treat=, covariates=) |
177.79 |
| backdoor, weighting | sp.ipw(...) |
-832.60, with a warning |
| front door | sp.front_door(..., outcome_model="additive") |
170.21 |
| instrument | sp.iv(...) |
first-stage F of 1.3: not usable |
The weighting estimate comes with a warning that the controls have an
effective sample size of 106 out of 963. The front door with the default
outcome_model="by_arm" raises, because nobody with low engagement won an
item and E[Y | D = 0, M = 1] has no data. The additive model extrapolates
under no interaction and should be reported as resting on that.
Refute¶
X = ["Guild Membership", "Player Skill Level", "Time Spent Playing"]
sp.refute(sp.aipw, df, y="In-game Purchases", treat="Side-quest Engagement",
covariates=X, method="placebo_treatment", seed=0).summary()
method= |
what is altered | the reruns should centre on |
|---|---|---|
"placebo_treatment" |
treatment shuffled across units | zero |
"dummy_outcome" |
outcome shuffled, or replaced by outcome_function(df) plus noise |
zero |
"random_common_cause" |
a random covariate added | the original estimate |
"data_subset" |
random subsets of the rows | the original estimate |
A shuffled outcome is a weak check, because shuffling removes the confounding along with the effect. The check that catches a missing adjustment keeps the confounding in place:
sp.refute(estimator, df, y=..., treat=..., method="dummy_outcome",
outcome_function=lambda f: 100 * f["Guild Membership"] + 50 * f["Player Skill Level"])
None of these detects an unobserved confounder. For that use
sp.sensemakr or sp.evalue.
Learn the graph¶
out = sp.pc_algorithm(df, ci_test="chi-square")
out["edges"], out["undirected_edges"], out["orientation_conflicts"]
The search is PC-stable and returns the skeleton and separating sets of
pcalg::pc. orientation_conflicts lists colliders the tests imply but
that could not be drawn because another collider had already oriented one
of their edges the other way. A non-empty list says the sample contradicts
itself about those edges.
Where StatsPAI differs from the book's tools¶
- The refutation p-value asks whether the reruns are centred where they should be, and a number of simulations too small to ever reject at the chosen level is refused.
sp.ipwreproduces a logit-based Horvitz-Thompson or Hajek estimate exactly. On the book's data that is not the figure DoWhy prints.- The front door fits each treatment arm separately by default. DoWhy's
two-stage version is
outcome_model="additive". - An implication that no stratum can test has a missing p-value, not a p-value of one.
The full account of the comparison is in
docs/dev/2026-10-06-ness-causal-ai-review.md.