Permutation invariance#
Give two people the same table with its columns in different orders. Ask each to learn a causal graph. They should get the same graph. Column order carries no information: it comes from a database export, a DataFrame, or whoever wrote the CSV. A method whose graph does not depend on that order is permutation invariant.
Many causal-discovery methods are not invariant. When two steps look equally good, they take the first one they find, and “first” means the leftmost column. Colombo and Maathuis (2014) showed this for PC. It matters for three reasons:
Reproducibility. Two analysts with the same data get different graphs.
Trust. An edge that flips when you reorder the columns was never supported by the data.
Benchmarks. Data generators often write columns in an order tied to the graph, so a method that favors the first column gains or loses accuracy it did not earn.
We plan to make Andrey the first package in which every method is permutation invariant. As a first step, 0.7.0 makes nine methods invariant: their choices follow the data, not the column order. Our benchmark also shuffles the columns, so their order cannot skew the results.
Andrey 0.7.0 brings permutation invariance to PC, FCI, GFCI, CD-NOD, hill climbing, tabu search, BOSS, GRaSP, and GIN: reordering the columns in your dataset no longer changes the discovered causal graph.
It costs almost no speed: each fix reuses work the method already does.
By eliminating sensitivity to arbitrary column ordering, Andrey makes causal discovery more reproducible, reliable, and trustworthy.
Datasets where reordering the columns changed the graph (lower is better)
| Andrey 0.7.0 | Andrey 0.6.0 | causal-learn | ||
|---|---|---|---|---|
| Score-based | Hill climbing, BDeu | 0% | 97% | - |
| Hill climbing, BIC | 0% | 98% | - | |
| Tabu search, BDeu | 0% | 75% | - | |
| Tabu search, BIC | 0% | 79% | - | |
| Permutation | BOSS | 0% | 97% | 91% |
| GRaSP | 0% | 99% | 100% | |
| Hidden causes | GIN | 0% | 73% | - |
How we found it#
It started with hill climbing. On the 0.6.0 categorical benchmark, pgmpy’s hill climbing had fewer
errors than Andrey’s on 17 datasets. We went looking for a bug and found a tie: adding x -> y and
adding y -> x can score exactly the same, and Andrey gave the tie to whichever column came first.
Reverse the columns, and Andrey had no more errors than pgmpy on 14 of the 17.
Then we checked our own data. Our generators wrote columns in an order tied to the graph: on a scale-free graph, the biggest hub came first. That order misled us. A better PC rule looked worse until we shuffled the columns.
Before 0.7.0: column 1 is the hub
0.7.0: columns shuffled
p-value rule minus column-order rule, endpoint errors
Columns as generated
Columns shuffled
fewer errorsmore errors
Why it happens#
Column order enters when a method takes the first of two equally good steps. In PC, two colliders can disagree about one edge, and whichever PC orients first wins.
Below, A causes B, C, and D, and each of them causes E. A and E look independent given B, with a
p-value of 0.052, just above the significance level of 0.05, so PC also proposes two false
colliders, A -> C <- E and A -> D <- E. The column-order rule visits colliders by the column
they point into, left to right: whenever C or D sits left of E, a false collider wins. The p-value
rule starts with the most clearly separated pair, so the true colliders always win.
Click a column to move it first, or shuffle them all.
The data. A causes B, C, and D, and each of them causes E. PC proposes five colliders; two are false.
0.6.0: column-order rule
Visits the colliders column by column, left to right: each waits under the column it points into. Click a column to move it first.
- 1A->C<-Ep = 0.052orients
- 2A->D<-Ep = 0.052orients
- 3B->E<-Cp = 0.821skipped:
C-Etaken - 4B->E<-Dp = 0.631skipped:
D-Etaken - 5C->E<-Dp = 0.273skipped:
C-E, D-Etaken
- 3B->E<-Cp = 0.821skipped:
C and D sit left of E, so their false colliders take C-E and D-E first.
This graph comes from 40 of 120 orders. Graphs seen: 1 of 4.
0.7.0: p-value rule
Visits the colliders by p-value, highest first: the most clearly separated pair goes first, wherever its columns sit.
- 1B->E<-Cp = 0.821orients
- 2B->E<-Dp = 0.631orients
- 3C->E<-Dp = 0.273orients
- 4A->C<-Ep = 0.052skipped:
C-Etaken - 5A->D<-Ep = 0.052skipped:
D-Etaken
The true graph, from all 120 orders.
The same thing happened in hill climbing and tabu search (score ties), BOSS and GRaSP (a starting order taken from the columns), GIN (a test built from each cluster’s first column), and FCI and GFCI (rules and a search step that visit the columns in order).
How we fixed it#
Each fix replaces a column’s position with a value computed from its data, which moves with the column.
PC. Conflicting colliders are oriented by the p-value that separated each pair, highest first, adapting Ramsey (2016). Meek’s rules then apply in rounds, so an edge they disagree about stays undirected. In 16 settings, no reordering changed PC’s graph, and PC made fewer errors in every one, for at most about 4% more time.
Hill climbing and tabu search. Exact score ties go to a rank computed from each column’s data. Accuracy is unchanged.
BOSS and GRaSP. Both search the columns sorted by that rank. Their accuracy matches causal-learn’s.
GIN. The test that orders the hidden causes uses every member of a cluster, not the first.
FCI and GFCI. Their rules apply in rounds, and the possible-d-separation step visits the columns by rank.
possible_dsep="fast"keeps the column-order step.Benchmark data. Our data generators and
andrey.data.load_datasetshuffle the columns along with the true graph.
Simulated data can understate the problem, so we also ran FCI on six real tables, each under 200 random column orders (20 for the full superconductivity table). With the default, no order changed the result.
Column orders that changed FCI's result (lower is better)
| Default | possible_dsep="fast" | ||
|---|---|---|---|
| Real tables | Sachs protein signaling11 variables, 853 rows | 0% | 0% |
| Wine13 variables, 178 rows | 0% | 0% | |
| Diabetes10 variables, 442 rows | 0% | 75% | |
| Breast cancer30 variables, 569 rows | 0% | 78% | |
| Superconductivity, 1,000-row sample81 variables, 1,000 rows | 0% | 64% | |
| Superconductivity81 variables, 21,263 rows | 0% | 100% | |
What it costs#
Almost nothing. A fix that runs more tests or more score evaluations would trade speed for order independence, so each fix reuses work the method already does:
PC orients by the p-values of the tests that removed each edge, so it runs no new test. With Meek’s rules in rounds, a fit takes at most about 4% longer.
Hill climbing and tabu search rank each column once per fit, from its data, and add no score evaluations.
BOSS and GRaSP sort the columns once by that rank and run the same search, in 1.03x and 0.91x the time of 0.6.0.
GIN runs its root test once per member of a cluster, but builds each column’s kernel matrix once per fit instead of once per test, so it runs about 4x to 7x faster than before.
FCI and GFCI put the columns in rank order once and run the same search, in 0.99x to 1.05x the time at 20 to 200 variables.
One trade remains. On some graphs with 400 or more variables, FCI’s and GFCI’s order-free step can
take longer: on one graph, 406 to 417 s against 120 to 336 s for the column-order step.
possible_dsep="fast" keeps that step.
How Andrey compares#
Permutation invariance is a yes-or-no property, so a few small examples answer it. We fit each package’s implementation at its defaults on three small datasets per method, each in nine column orders, and checked whether any reordering changed the graph.
Does reordering the columns change the graph?
| Andrey 0.7.0 | causal-learn 0.1.4.8 | pcalg 2.7-12 | bnlearn 5.2.1 | pgmpy 1.1.2 | Tetrad 7.6.10 | ||
|---|---|---|---|---|---|---|---|
| Constraint | PC | Same | Changed | Changed | Changed | By name | By name |
| FCI | Same | Changed | Changed | - | - | By name | |
| Score-based | Hill climbing | Same | - | - | Changed | By run | - |
| Tabu search | Same | - | - | Changed | - | - | |
| Permutation | BOSS | Same | Changed | - | - | - | Changed |
| GRaSP | Same | Changed | - | - | - | Changed | |
| Hidden causes | GIN | Same | Changed* | - | - | - | Changed |
For BOSS, GRaSP, GIN, and tabu search, every other implementation we checked changed its graph under some reordering. pcalg and causal-learn can make PC’s colliders order-independent, but only with options that are off by default; Andrey does it by default, with no extra tests.
Try it#
pip install -U andrey-core
Fit PC on 10 random orders of the same columns, with the old rule and with 0.7.0’s default:
import numpy as np
import andrey
from andrey.data import sample_scm
data = sample_scm(d=50, n=500, density=4.0, scale="standardize", seed=0).data
rng = np.random.default_rng(0)
orders = [rng.permutation(data.shape[1]) for _ in range(10)]
for rule in ("column", "p_value"):
graphs = set()
for order in orders:
fit = andrey.pc(data[:, order], alpha=0.05, collider_order=rule)
back = np.argsort(order) # put the columns back in their first order
graphs.add(fit.structure.adjacency[np.ix_(back, back)].tobytes())
count = len(graphs)
print(f"collider_order={rule!r}: {count} graph{'s' if count > 1 else ''} from 10 column orders")
collider_order='column': 10 graphs from 10 column orders
collider_order='p_value': 1 graph from 10 column orders
The 0.7.0 release card reports every measurement in this release, with each method’s limits, and the changelog lists every change.
Andrey is in alpha. The API may change between releases. When you try it on your own data, tell us what breaks or which method you need. We also welcome pull requests. The contributing guide helps you get set up.