GuidesDecode a DAG from edge probabilities

Decode a DAG from edge probabilities

Turn a matrix of edge probabilities into a DAG with a named decoder, and see which edges it dropped and why. You'll compare decoders, draw candidate DAGs with a seeded sampler, and constrain the result with domain knowledge.

Run the program#

The matrix is the one from Work with edge probabilities: four edges pass 0.5, and three of them form a cycle.

import numpy as np

from ergodic import DomainKnowledge, EdgeProbabilities, decode

# 1. Independent edge probabilities; three edges above 0.5 form a cycle.
names = ("rain", "sprinkler", "wet", "slippery")
values = np.array(
    [
        [0.00, 0.45, 0.91, 0.06],
        [0.40, 0.00, 0.83, 0.05],
        [0.04, 0.07, 0.00, 0.88],
        [0.58, 0.02, 0.10, 0.00],
    ]
)
probs = EdgeProbabilities(names, values, kind="independent")

# 2. The default decoder: greedy pair-MAP.
result = decode(probs)
print(result)
print(result.graph)
print(result.diagnostics["rejected_cycles"], result.diagnostics["note"])

# 3. A lower threshold nominates more edges.
print(decode(probs, threshold=0.3).graph)

# 4. Exact MAP and the expected-count decoder.
print(decode(probs, "exact_map").graph)
print(decode(probs, "expected_count").graph)

# 5. Candidate DAGs from a seeded sampler, and their consensus.
sampled = decode(probs, "pair_sampler", samples=200, seed=0)
print(sampled)
print(sampled.graph)

# 6. Knowledge: nothing causes the sprinkler.
knowledge = DomainKnowledge().root("sprinkler")
print(decode(probs, threshold=0.3, knowledge=knowledge).graph)

DecodeResult(decode:greedy_v1: decoded dag, 4 nodes, 3 edges)
MixedGraph(kind=dag, nodes=4, edges=[rain --> wet, slippery <-- wet, sprinkler --> wet])
(1,) not a graph posterior; orientation not certified
MixedGraph(kind=dag, nodes=4, edges=[rain --> sprinkler, rain --> wet, slippery <-- wet, sprinkler --> wet])
MixedGraph(kind=dag, nodes=4, edges=[rain --> wet, slippery <-- wet, sprinkler --> wet])
MixedGraph(kind=dag, nodes=4, edges=[rain --> sprinkler, rain --> wet, slippery <-- wet, sprinkler --> wet])
DecodeResult(decode:pair_sampler_v1: decoded dag, 4 nodes, 3 edges, consensus of 200 candidates)
MixedGraph(kind=dag, nodes=4, edges=[rain --> wet, slippery <-- wet, sprinkler --> wet])
MixedGraph(kind=dag, nodes=4, edges=[rain <-- sprinkler, rain --> wet, slippery <-- wet, sprinkler --> wet])

• Step 2. Greedy nominates each pair's likelier direction when it beats "no edge", which on independent marginals means above 0.5. It adds them strongest first and skips any that would close a cycle: slippery --> rain, the weakest of the three, is the one rejected cycle. The result is a DecodeResult, a DiscoveryResult labelled as decoded.
• Step 3. With threshold=0.3 an edge is nominated when its own value is above 0.3, so rain --> sprinkler (0.45) joins. The cut is strict, as AVICI and CDFM apply theirs.
• Step 4. Exact MAP finds the DAG with the best total score and agrees with greedy here. The expected-count decoder keeps as many edges as the probabilities expect, rounded (four), so it adds rain --> sprinkler too.
• Step 5. The sampler draws 200 DAGs. graph is their medoid, the candidate closest to all the others (the earliest one on an exact tie); every candidate is in candidates.
• Step 6. root("sprinkler") forbids every edge into sprinkler, so the pair flips to sprinkler --> rain, its other direction at 0.40.

Choose a decoder#

Every decoder reads the pair-state view (see Work with edge probabilities) and scores an edge by w(i -> j) = log P(i -> j) - log P(no edge).

Greedy is the default. In the research it came from it matched exact MAP on 238 of 240 graphs. It isn't best everywhere: on one physical benchmark the expected-count decoder scored better. The decoders that search orderings exactly refuse more than max_nodes variables (20 by default), because that search can grow exponentially.

Read the result#

decode returns a DecodeResult. graph is a DAG over every variable, method is "decode:<id>", candidates holds every decoded DAG and probabilities the input. diagnostics records:

Decoding is deterministic: ties go to the earlier variable, the same seed gives the same candidates, and neither the input nor numpy's global random state changes.

Use domain knowledge#

decode(..., knowledge=...) applies hard constraints:

• Forbidden edges, tiers, roots, sinks and restrict_to remove directed edges before decoding.
• Order (order, tiers) is enforced as one time order: the decoded DAG has a single topological order that agrees with every order statement and every edge at once. For tiers that's the same as "no later variable is an ancestor of an earlier one". For separate order() statements it's stricter than DomainKnowledge.violations, which checks each statement on its own: with order("a", "b") and order("c", "d"), the edges b --> c and d --> a pass violations together, but no single order fits them, so a decoder keeps only one.
• Required edges go in first, whatever their probability. A required set with a cycle raises.
• confound raises: a DAG can't hold a latent confounder.

Every variable the knowledge names must be in the probabilities. The decoded graph is checked against the knowledge before it's returned.

Decode a Welt model's scores#

AVICI and CDFM return directed graphs that may hold cycles, so ergodic.welt.discover refuses them before anything is paid for. Pass decoder= to opt in: discover(data, "cdfm", decoder="greedy") runs the model once, reads its scores as independent edge probabilities and decodes them. Greedy then uses the model's own threshold (0.5 for AVICI, CDFM's adaptive one). The result keeps the scores in diagnostics["edge_probabilities"], so you can try another decoder locally without a new run. The exact decoders' max_nodes cap (20 by default) is checked against your data's variables before anything is uploaded. See Decode AVICI and CDFM.

Arrow returns a DAG of its own, so ergodic uses that by default. Its scores are decoded only if you ask, with decoder=, and are then read as independent marginals, like AVICI's and CDFM's. That's a modelling choice: Welt doesn't state that Arrow's scores are marginals.

Next steps#