Causal Discovery

Python · causal-learn  ·  R · bnlearn  ·  simulation with a known DAG, and Sachs et al. (2005)

What Observational Data Contains

The previous page assumed the graph was known and read identification off it. That raises the obvious question: where does the graph come from? Domain knowledge, usually. But the structure itself can be learned from data, and what that learning can and cannot deliver is unusually well characterised.

The central limit is Markov equivalence. Several distinct DAGs can imply exactly the same set of conditional independencies, and observational data sees nothing but those independencies — so it cannot tell those graphs apart. Discovery algorithms therefore return not a DAG but a CPDAG: some edges directed, others left undirected because the data does not contain the answer.

The PC algorithm works by constraint: test conditional independencies to recover the skeleton — which variables are directly connected, before any arrow has a direction — then orient v-structures — the unshielded colliders, where two arrows meet at a common effect that connects them to nothing else — and propagate whatever follows logically from those. Simulated from a known five-variable DAG with Gaussian noise, it recovers the entire skeleton, orients that collider X0X2X1X_0 \rightarrow X_2 \leftarrow X_1 and the downstream X2X3X_2 \rightarrow X_3 — and leaves X0X4X_0 - X_4 undirected. That edge is not a failure of the algorithm. Both orientations imply identical independencies, so no amount of Gaussian observational data distinguishes them.

GES takes the other route — greedy search over graphs scoring by BIC rather than testing independencies — and returns the same CPDAG, including the same undirected edge. Two different algorithm families reaching the same boundary is the point: it is a property of the data, not of either method.

Non-Gaussianity Breaks the Tie

LiNGAM escapes it with one extra ingredient. If the noise is non-Gaussian, an asymmetry appears that Gaussianity hides: for Y=βX+εY = \beta X + \varepsilon, the regression residual is independent of the regressor only in the true causal direction. Re-running the same structure with uniform noise, LiNGAM recovers every direction, orienting X4X0X_4 \rightarrow X_0 correctly — the edge PC could not touch.

The edge list is read off the fitted coefficient matrix with a 0.15 cutoff, which looks like a tuning knob and is not one. DirectLiNGAM prunes to exact zeros: the fitted matrix has four non-zero entries, the smallest 0.63 and the largest 0.91, with nothing in between and nothing below. Any cutoff from 0.05 to 0.40 returns the identical graph, and the recovered coefficients — 0.63, 0.70, 0.79, 0.91 — sit on the true 0.6, 0.7, 0.8, 0.9.

All of it rests on assumptions that cannot be checked from the data being analysed. Causal sufficiency — no unmeasured common causes — is the load-bearing one, and its failure is not subtle: two variables with a hidden common cause and no edge between them produce a phantom edge, which PC duly reports. Faithfulness requires that the only independencies in the data are those the graph implies, ruling out exact cancellations. And the independence tests must be reliable in finite samples, which is where the real data comes in.

The Benchmark, and What It Actually Returns

Sachs et al. (2005) measured 11 phosphorylated proteins in thousands of individual human immune cells by flow cytometry. It is the field's benchmark because decades of molecular biology supply an established ground-truth network — the Raf→Mek→Erk cascade, PKA and PKC as broad regulators, the Plcg→PIP3→PIP2 messenger chain — against which a recovered graph can actually be scored. Running PC on the 853 observational cells and comparing skeletons:

PC against the consensus networkresult
true edges recovered7 of 17
false edges0
precision1.00
recall0.41

853 observational cells, 11 proteins, skeleton comparison. The consensus network has 17 directed edges.

Precision is not high here. It is perfect — not one false edge — and that deserves an explanation rather than an adjective, because the explanation is unflattering. PC asserted 7 adjacencies out of 55 possible pairs, on data whose largest off-diagonal correlation is 0.82. It is not being accurate so much as being reluctant, and what little it commits to is safe.

Which raises the question the write-up usually skips: are the missing edges near misses? The test level was fixed at 0.05 without comment. Sweeping it:

αedges foundtruefalseprecisionrecall
0.0016601.000.35
0.016601.000.35
0.05 the default7701.000.41
0.107701.000.41
0.2010730.700.41

Twenty times more permissive, three more false edges, one more true one.

The recall column does not move. Loosening the test twentyfold, from 0.01 to 0.20, adds four edges of which three are false, and takes true positives from 6 to 7. The ten missing edges are not sitting just past a significance boundary waiting for a more permissive threshold — they are invisible to conditional-independence testing on this data, and no choice of α\alpha recovers them. That is a statement about what the observational distribution contains, not about how the test was tuned.

The R companion is a genuine cross-check rather than a restatement. bnlearn's pc.stable and its score-based hill-climbing both return the same seven edges — TP 7, FP 0, recall 0.41 — on a different implementation, in a different language, with a different algorithm family for the second. The boundary is in the data.

One preprocessing note, since the log transform was taken as standard rather than checked. On raw abundances PC finds 8 edges, still with no false ones, for a recall of 0.47 against the logged 0.41. The transform is conventional for flow cytometry and it costs a true edge here. Small, and worth knowing that the conventional choice was a choice.

The honest stance is that discovery is a hypothesis generator, not a design. It narrows the space of graphs worth considering, which is genuinely useful when the space is large and domain knowledge is thin. It does not hand you a network, and on the field's own benchmark it recovers under half of one. The ending is the arc's through-line stated by a different discipline: Sachs et al. reconstructed the full directed network by intervening — perturbing each protein in turn — not by analysing the observational data harder. The canonical demonstration that structure can be learned from data was settled by running experiments.

Where this sits

This page is the inverse of DAGs & the SCM: there the graph was given and identification read off it, here only its equivalence class is recoverable. Markov equivalence is the structural face of d-separation, and causal sufficiency is the discovery version of the unconfoundedness assumption — untestable in exactly the same way, and failing into a phantom edge rather than a biased coefficient. That the benchmark was completed by intervention rather than observation echoes the randomized-experiments group, where the arc began.

Notebooks

Downloads

References