Analysis

jaxonomy.analysis

Whole-model analysis on top of the framework's dependency structure.

The influence graph merges the model's leaf-level dependency DAG (which says whether information flows) with autodiff Jacobians (which say how much), giving quantitative model slicing, chain-rule path attribution, bottleneck detection, and dead-edge diagnostics on one queryable object.

InfluenceGraph dataclass

A model's dependency structure with autodiff-computed edge weights.

Build with :func:influence_graph; see this module's docstring for the weighting conventions the numbers obey.

Attributes:

Name	Type	Description
system	Any	The analyzed Diagram or LeafSystem.
graph	DiGraph	The underlying networkx.DiGraph. Node attributes describe the signal (kind, block, port, size, value, units, sample_time, hybrid); edge attributes carry kind, jacobian, relative, weight, magnitude, local_gradient, note and — in trajectory mode — profile.
tau	float	Time scale applied to continuous-state-rate edges, in seconds.
normalize	str	"relative" or "none".
scale_floor	float	Lower bound on a signal's magnitude when normalizing.
at	str	"operating_point" or "trajectory".
times	Optional[ndarray]	Snapshot times in trajectory mode, else None.
reduce	str	How a trajectory profile became the scalar weight.
block_notes	Dict[str, Dict[str, str]]	Per-block explanations for anything not differentiated.

Source code in jaxonomy/analysis/influence.py

@dataclass
class InfluenceGraph:
    """A model's dependency structure with autodiff-computed edge weights.

    Build with :func:`influence_graph`; see this module's docstring for the
    weighting conventions the numbers obey.

    Attributes:
        system: The analyzed ``Diagram`` or ``LeafSystem``.
        graph: The underlying ``networkx.DiGraph``. Node attributes describe
            the signal (``kind``, ``block``, ``port``, ``size``, ``value``,
            ``units``, ``sample_time``, ``hybrid``); edge attributes carry
            ``kind``, ``jacobian``, ``relative``, ``weight``, ``magnitude``,
            ``local_gradient``, ``note`` and — in trajectory mode —
            ``profile``.
        tau: Time scale applied to continuous-state-rate edges, in seconds.
        normalize: ``"relative"`` or ``"none"``.
        scale_floor: Lower bound on a signal's magnitude when normalizing.
        at: ``"operating_point"`` or ``"trajectory"``.
        times: Snapshot times in trajectory mode, else None.
        reduce: How a trajectory profile became the scalar weight.
        block_notes: Per-block explanations for anything not differentiated.
    """

    system: Any
    graph: nx.DiGraph
    tau: float
    normalize: str
    scale_floor: float
    at: str
    times: Optional[np.ndarray]
    reduce: str
    block_notes: Dict[str, Dict[str, str]] = field(default_factory=dict)
    structure: Optional[nx.DiGraph] = None

    # -- basics ------------------------------------------------------------

    def __repr__(self) -> str:
        return (
            f"InfluenceGraph({self.system.name}, {self.graph.number_of_nodes()} nodes, "
            f"{self.graph.number_of_edges()} edges, at={self.at!r}, tau={self.tau:g})"
        )

    @property
    def n_blocks(self) -> int:
        return len({d["block"] for _, d in self.graph.nodes(data=True)})

    @property
    def blocks(self) -> List[str]:
        return sorted({d["block"] for _, d in self.graph.nodes(data=True)})

    def resolve(self, spec) -> str:
        """Turn a port object, locator, or name fragment into a node id.

        Accepts an exact node id, an ``InputPort`` / ``OutputPort``, a
        ``(system, port_index)`` locator, or any unambiguous suffix of a node
        id (``"integ:out:out_0"``, ``"integ"``, ``"out:y"``).
        """
        if isinstance(spec, (InputPort, OutputPort)):
            return self._require(port_node_id(spec))
        if isinstance(spec, tuple) and len(spec) == 2:
            system, index = spec
            ports = system.output_ports if hasattr(system, "output_ports") else []
            if index < len(ports):
                return self._require(port_node_id(ports[index]))
        if not isinstance(spec, str):
            raise TypeError(
                f"Cannot resolve {spec!r} to an influence-graph node; pass a node "
                f"id, a port object, or a (system, port_index) locator."
            )
        if spec in self.graph:
            return spec
        matches = [n for n in self.graph if n.endswith(spec) or spec in n]
        if len(matches) == 1:
            return matches[0]
        if not matches:
            raise KeyError(
                f"No influence-graph node matches {spec!r}. "
                f"Nodes are named '<block path>:in|out:<port>' or "
                f"'<block path>:xc|xd'; {self.graph.number_of_nodes()} exist, "
                f"e.g. {sorted(self.graph)[:4]}."
            )
        raise KeyError(
            f"{spec!r} is ambiguous — it matches {len(matches)} nodes: "
            f"{sorted(matches)[:6]}. Use a longer fragment or the full node id."
        )

    def _require(self, node_id: str) -> str:
        if node_id not in self.graph:
            raise KeyError(
                f"{node_id!r} is not in this influence graph. It may belong to a "
                f"different model, or be a port of a sub-Diagram rather than a leaf."
            )
        return node_id

    # -- traversal ---------------------------------------------------------

    def _step_weight(self, src: str, dst: str) -> Tuple[float, bool]:
        """(magnitude used for propagation, whether it is actually known)."""
        data = self.graph.edges[src, dst]
        if not data["local_gradient"]:
            # Unknown, not zero. The unit gain is a placeholder that keeps a
            # product finite for reporting; callers must branch on the second
            # element rather than treat this as a measurement.
            return 1.0, False
        return abs(data["magnitude"]), True

    def _amplification(
        self, direction: str, max_depth: int
    ) -> Tuple[List[Dict[str, float]], List[Dict[str, bool]]]:
        """The pruning bound: how much a continuation could still contribute.

        Returns ``(bounds, unmeasurable)``, both indexed by remaining depth.
        ``bounds[d][n]`` is the largest product over *measurable* paths of ≤ ``d``
        edges arriving at ``n``; ``unmeasurable[d][n]`` says whether any such
        path crosses an edge with no local gradient.

        The bound exists because a running product is not an upper bound on the
        finished path's — a relative weight is an elasticity and can exceed 1, so
        a partial path can dip below any threshold and climb back. Relaxing over
        depth gives the largest factor still available, in ``max_depth`` sweeps
        of the edge list, and being an over-estimate it can only prune paths that
        could not have qualified.

        The two results are kept separate on purpose. An unmeasurable
        continuation must never be pruned — "we could not compute it" is not
        "it is small" — but folding that into the *numeric* bound as an infinite
        gain would propagate infinity to every node upstream of any comparator
        or quantizer, switching pruning off across the whole model and
        collapsing the search back to an exponential one. A boolean reachability
        flag suppresses pruning exactly where an unknown edge is actually
        within reach, and leaves the numeric bound finite everywhere else.
        """
        incoming = (
            self.graph.predecessors if direction == "backward" else self.graph.successors
        )
        bounds = [{node: 1.0 for node in self.graph}]
        unmeasurable = [{node: False for node in self.graph}]
        for _ in range(max_depth):
            previous_bound, previous_flag = bounds[-1], unmeasurable[-1]
            current_bound, current_flag = dict(previous_bound), dict(previous_flag)
            for node in self.graph:
                for other in incoming(node):
                    if other == node:
                        continue
                    edge = (
                        (other, node) if direction == "backward" else (node, other)
                    )
                    magnitude, known = self._step_weight(*edge)
                    if not known:
                        current_flag[node] = True
                        continue
                    candidate = previous_bound[other] * magnitude
                    if candidate > current_bound[node]:
                        current_bound[node] = candidate
                    if previous_flag[other]:
                        current_flag[node] = True
            if current_bound == previous_bound and current_flag == previous_flag:
                break
            bounds.append(current_bound)
            unmeasurable.append(current_flag)
        return bounds, unmeasurable

    def _reach(
        self,
        target: str,
        max_depth: int,
        direction: str,
        threshold: float = 0.0,
        max_expansions: int = 200_000,
        amplification: Optional[
            Tuple[List[Dict[str, float]], List[Dict[str, bool]]]
        ] = None,
    ) -> Tuple[
        Dict[str, float], Dict[str, bool], Dict[str, Tuple[str, ...]], bool
    ]:
        """Best *simple*-path product between every node and ``target``.

        Three properties this has to get right, each of which a more obvious
        implementation gets wrong on a real model:

        **Simple paths only.** Allowing a path to revisit a node lets it
        circulate a feedback loop, multiplying by the loop gain each turn, and a
        node's "influence" then reports how many times the search went round
        rather than how much signal gets through — in a closed loop it even
        makes a node influence *itself* by a large factor. Restricting to simple
        paths is also the standard signal-flow-graph notion of a forward path,
        and matches what :meth:`attribute` enumerates.

        **No pruning on the running product.** A relative weight is an
        elasticity, not a gain bounded by 1: a summing junction whose output
        nearly cancels (any controller error signal near steady state) has an
        elasticity far above 1, so a partial product can dip below any threshold
        and climb back above it. Pruning uses :meth:`_amplification` instead,
        which bounds what a continuation could still contribute.

        **Unmeasurable edges are reachability, not search.** Past an edge with
        no local gradient there is no product left to maximize, and *before* one
        there is nothing to optimize either — only the question of what connects.
        Both are answered by two linear BFS sweeps after the path search, so a
        comparator or quantizer anywhere upstream costs a pass over the edge
        list rather than an unpruned walk of everything between it and the
        target.

        Returns ``(best, unknown, routes, truncated)``, where ``routes[n]`` is
        the node sequence realizing ``best[n]`` (from ``n`` to ``target``, or
        the reverse for a forward search) and ``truncated`` is True when the
        expansion budget ran out, in which case the scores are lower bounds.
        Scores for nodes the sweeps supplied are lower bounds by construction —
        one real route's product rather than the best one's.
        """
        step = (
            self.graph.predecessors if direction == "backward" else self.graph.successors
        )
        if amplification is None and threshold > 0.0:
            amplification = self._amplification(direction, max_depth)
        best: Dict[str, float] = {target: 1.0}
        unknown: Dict[str, bool] = {target: False}
        # The route realizing each node's best product. Slicing needs it to stay
        # connected: retaining nodes by their own influence alone leaves the
        # intermediate hops of a dominant route out, and a disconnected slice
        # makes `subgraph` and `bottlenecks` meaningless.
        routes: Dict[str, Tuple[str, ...]] = {target: (target,)}
        truncated = False
        expansions = 0

        # Path enumeration runs over fully measurable routes only, and prunes
        # on the numeric bound alone. Everything to do with unmeasurable edges
        # is handled by the two reachability sweeps below, because *both* halves
        # of that question — getting to an unknown edge and going past it — are
        # reachability, not optimization. Letting the DFS off its leash to find
        # the unknown frontier (the obvious alternative) means walking a dense
        # region as an unpruned simple-path enumeration before ever arriving at
        # the comparator that made it necessary.
        stack = [(target, 1.0, frozenset((target,)), 0, (target,))]
        while stack:
            node, product, trail, depth, route = stack.pop()
            if depth >= max_depth:
                continue
            expansions += 1
            if expansions > max_expansions:
                truncated = True
                break
            for other in step(node):
                if other in trail:
                    continue  # a repeat would be a loop turn, not a new path
                edge = (other, node) if direction == "backward" else (node, other)
                magnitude, known = self._step_weight(*edge)
                if not known:
                    continue  # picked up by the frontier sweep below
                new_product = product * magnitude
                new_route = (other,) + route
                if amplification is not None:
                    bounds, _ = amplification
                    remaining = min(max_depth - depth - 1, len(bounds) - 1)
                    if new_product * bounds[remaining][other] < threshold:
                        continue
                previous = best.get(other)
                if previous is None or new_product > previous:
                    best[other] = new_product
                    routes[other] = new_route
                stack.append(
                    (other, new_product, trail | {other}, depth + 1, new_route)
                )

        # Sweep 1: how far is each node from the target, by what route, and what
        # that route carries. Plain BFS, so it costs one pass whatever the
        # model's density. The product is needed because a seed's route runs
        # through nodes the DFS pruned; without a score of their own they would
        # be dropped from the slice while their edges were kept, leaving the
        # retained set inconsistent with the subgraph.
        hops: Dict[str, int] = {target: 0}
        to_target: Dict[str, Tuple[str, ...]] = {target: (target,)}
        hop_product: Dict[str, float] = {target: 1.0}
        hop_unknown: Dict[str, bool] = {target: False}
        frontier_queue = deque([target])
        while frontier_queue:
            node = frontier_queue.popleft()
            if hops[node] >= max_depth:
                continue
            for other in step(node):
                if other == node or other in hops:
                    continue
                edge = (other, node) if direction == "backward" else (node, other)
                magnitude, known = self._step_weight(*edge)
                hops[other] = hops[node] + 1
                to_target[other] = (other,) + to_target[node]
                hop_product[other] = hop_product[node] * magnitude
                hop_unknown[other] = hop_unknown[node] or not known
                frontier_queue.append(other)

        # Sweep 2: seed at every unmeasurable edge that can still reach the
        # target, then walk away from it. This is what keeps a comparator, a
        # quantizer, or anything behind them in the answer.
        unknown_seeds: List[Tuple[str, int, float, Tuple[str, ...]]] = []
        for source_node, destination, data in self.graph.edges(data=True):
            if data["local_gradient"] or source_node == destination:
                continue
            head, tail = (
                (destination, source_node)
                if direction == "backward"
                else (source_node, destination)
            )
            reached = hops.get(head)
            if reached is None or reached + 1 > max_depth:
                continue
            # Everything the seed's route passes through has to carry a score
            # too, or it would be filtered out of the slice while its edges
            # stayed in. The BFS route is one real path, so its product is a
            # lower bound on that node's influence — honest, and finite.
            for hop in to_target[head]:
                best.setdefault(hop, hop_product[hop])
                routes.setdefault(hop, to_target[hop])
                if hop_unknown[hop]:
                    unknown[hop] = True
            route = (tail,) + to_target[head]
            unknown[tail] = True
            best.setdefault(tail, 1.0)
            routes.setdefault(tail, route)
            unknown_seeds.append((tail, reached + 1, best[tail], route))

        # Scores past an unmeasurable edge are placeholders, not measurements —
        # `unknown` marks them so no caller reads them as one.
        queue = deque(unknown_seeds)
        shallowest: Dict[str, int] = {}
        while queue:
            node, depth, product, route = queue.popleft()
            if depth >= max_depth or shallowest.get(node, max_depth + 1) <= depth:
                continue
            shallowest[node] = depth
            for other in step(node):
                if other == node or other in route:
                    continue
                unknown[other] = True
                best.setdefault(other, product)
                extended = (other,) + route
                routes.setdefault(other, extended)
                queue.append((other, depth + 1, product, extended))

        return best, unknown, routes, truncated

    def slice(
        self,
        target,
        threshold: float = 0.01,
        *,
        direction: str = "backward",
        max_depth: int = 32,
    ) -> InfluenceSlice:
        """Quantitative model slice: what influences ``target`` by ≥ ``threshold``.

        The boolean answer — everything structurally upstream — is
        :meth:`structural_slice`; this one keeps only what lies on a path
        carrying at least ``threshold`` of the influence, in the
        relative-sensitivity sense described in the module docstring.

        ``0.01`` reads as "1%" only when ``tau`` is comparable to the time
        constants on the paths involved — the threshold is absolute, and a
        path across *k* integrators carries a factor of ``tau**k``, so the
        same cutoff means different things at different ``tau``. When the
        strongest contributor scores 95, ``threshold=0.01`` retains everything
        down to ~0.01% of it, not 1%. Use :meth:`relative_threshold` to get
        the cutoff that means a fraction *of the dominant contributor*.

        Two kinds of node are kept, and the distinction is load-bearing. A node
        is **influential** when its own best path to ``target`` clears the
        threshold. It is a **connector** when it merely lies on some influential
        node's best route: a relative weight is an elasticity, so a signal can
        pass through a junction that nearly cancels it and be amplified back
        afterwards, leaving a mid-route node with a small score of its own.
        Keeping only the influential ones would punch holes in the result —
        naming a block as influential while the route from it to the target ran
        through blocks that had been dropped, leaving
        :attr:`InfluenceSlice.subgraph` disconnected and :meth:`bottlenecks`
        meaningless. Connectors are read off the routes the search actually
        found, so nothing is added that no real path uses.

        ``scores`` reports every retained node's own best product to the target,
        which is the number to rank by.

        Args:
            target: Node id, port object, or name fragment (see :meth:`resolve`).
            threshold: Minimum ``|path product|`` for a path to be retained.
            direction: ``"backward"`` (default, what influences the target) or
                ``"forward"`` (what the target influences).
            max_depth: Hard bound on path length, and the only hard bound — a
                partial product is not a bound on the whole path's (see
                :meth:`_reach`), so nothing may be pruned on the running value.

        Returns:
            An :class:`InfluenceSlice`.
        """
        if direction not in ("backward", "forward"):
            raise ValueError(
                f"direction must be 'backward' or 'forward', got {direction!r}"
            )
        node = self.resolve(target)
        best, unknown, routes, truncated = self._reach(
            node, max_depth, direction, threshold=threshold
        )

        influential = {
            other
            for other, score in best.items()
            if other == node or score >= threshold or unknown.get(other, False)
        }
        # Pull in whatever each influential node's own best route passes
        # through, so the result is a connected sub-model rather than a set of
        # names with no way to get between them.
        retained_nodes = set(influential)
        for other in influential:
            retained_nodes.update(routes.get(other, (other,)))

        # Every node on a retained route gets a score. Dropping the ones the
        # numeric search never scored would leave `blocks` naming fewer blocks
        # than `edges` actually connects — the retained set and the subgraph
        # would disagree, and the route's interior blocks would vanish from the
        # answer even though they are the only way the influence travels.
        retained = {
            other: best[other] if other in best else 0.0 for other in retained_nodes
        }
        retained.setdefault(node, best[node])

        kept = set()
        for other in influential:
            route = routes.get(other, ())
            for first, second in zip(route, route[1:]):
                kept.add(
                    (first, second) if direction == "backward" else (second, first)
                )
        # Beyond the best routes, keep any edge between retained nodes that
        # itself carries a qualifying path — a parallel branch of comparable
        # strength belongs in the slice even though some other route was best.
        for other in retained:
            neighbours = (
                self.graph.successors(other)
                if direction == "backward"
                else self.graph.predecessors(other)
            )
            for downstream in neighbours:
                if downstream == other or downstream not in retained:
                    continue
                edge = (
                    (other, downstream)
                    if direction == "backward"
                    else (downstream, other)
                )
                magnitude, known = self._step_weight(*edge)
                if not known or magnitude * retained[downstream] >= threshold:
                    kept.add(edge)

        blocks = sorted({self.graph.nodes[n]["block"] for n in retained})
        return InfluenceSlice(
            target=node,
            threshold=threshold,
            direction=direction,
            scores=dict(retained),
            edges=sorted(kept),
            blocks=blocks,
            unknown_nodes=sorted(
                other for other in retained if unknown.get(other, False)
            ),
            graph=self,
            truncated=truncated,
        )

    def relative_threshold(
        self,
        target,
        fraction: float = 0.01,
        *,
        direction: str = "backward",
        max_depth: int = 32,
        floor: float = 1e-12,
    ) -> float:
        """A threshold set at ``fraction`` of the strongest influence on ``target``.

        An absolute threshold only reads as a percentage when ``tau`` is
        comparable to the time constants on the paths involved (see the module
        docstring). Scaling to the strongest score makes "keep what carries at
        least 1% of what the dominant contributor carries" mean the same thing at
        any ``tau``.

        Args:
            target: Node id, port object, or name fragment.
            fraction: Fraction of the strongest score to keep.
            direction: As in :meth:`slice`.
            max_depth: As in :meth:`slice`.
            floor: Threshold the reference sweep runs at, and the value returned
                when nothing upstream carries influence. It is passed to the
                search rather than left at zero so the sweep stays pruned; a
                model whose strongest contributor falls below it would yield
                ``floor`` itself.

        Returns:
            A threshold to pass to :meth:`slice` / :meth:`bottlenecks`.
        """
        node = self.resolve(target)
        best, _, _, _ = self._reach(node, max_depth, direction, threshold=floor)
        others = [value for other, value in best.items() if other != node]
        if not others or max(others) <= 0.0:
            return floor
        return max(fraction * max(others), floor)

    def structural_slice(self, target, *, direction: str = "backward") -> List[str]:
        """Boolean slice: every block structurally connected to ``target``.

        The over-approximation :meth:`slice` improves on, computed from the
        model's declared connectivity rather than from the weighted graph — so
        it stays a genuine bound even where a Jacobian could not be taken.
        Provided so the two can be compared directly on a real model.
        """
        node = self.resolve(target)
        structure = self.structure if self.structure is not None else self.graph
        if node not in structure:
            return [self.graph.nodes[node]["block"]]
        reached = (
            nx.ancestors(structure, node)
            if direction == "backward"
            else nx.descendants(structure, node)
        )
        return sorted({structure.nodes[n]["block"] for n in reached | {node}})

    def _enumerate_paths(
        self,
        source: str,
        target: str,
        threshold: float,
        max_depth: int,
        max_paths: int,
    ) -> Tuple[List[Dict[str, Any]], bool]:
        paths: List[Dict[str, Any]] = []
        truncated = False

        # The running product alone is not a bound on the finished path's — an
        # elasticity above 1 further along can lift it back over the threshold
        # (see :meth:`_reach`). Pruning is therefore done against
        # running x (the most any continuation could still contribute), which is
        # admissible and so never discards a path that would have qualified.
        #
        # This bound comes from :meth:`_amplification`, not from :meth:`_reach`:
        # `_reach` is itself a simple-path DFS, and calling it here without a
        # threshold would be an unpruned exponential search whose expansion
        # budget could silently cut off nodes — making this enumeration report
        # "no path" while claiming it was complete.
        onward_bounds, onward_unmeasurable = self._amplification(
            "forward", max_depth
        )

        # Pruning is switched off once a path turns unknown (its product is a
        # placeholder, so the bound says nothing), which leaves the walk through
        # an unmeasurable region unbounded. `max_paths` only fires when paths are
        # actually found, so a region that reaches the target rarely needs this
        # second cap.
        expansions = 0
        max_expansions = 200_000

        stack = [(source, [source], 1.0, True, False)]
        while stack:
            node, trail, product, signed, unknown = stack.pop()
            if node == target and len(trail) > 1:
                paths.append(
                    {
                        "nodes": list(trail),
                        "product": product,
                        "signed": signed,
                        "unknown": unknown,
                    }
                )
                if len(paths) >= max_paths:
                    truncated = True
                    break
                continue
            if len(trail) > max_depth:
                truncated = True
                continue
            expansions += 1
            if expansions > max_expansions:
                truncated = True
                break
            for successor in self.graph.successors(node):
                if successor in trail:
                    continue  # a loop turn cannot add a new path
                data = self.graph.edges[node, successor]
                if not data["local_gradient"]:
                    stack.append(
                        (successor, trail + [successor], product, False, True)
                    )
                    continue
                weight = data["weight"]
                block = data["relative"]
                step_signed = signed and block is not None and block.size == 1
                new_product = product * (weight if step_signed else abs(weight))
                if not unknown:
                    depth_left = min(max_depth - len(trail), len(onward_bounds) - 1)
                    if depth_left < 0:
                        continue
                    # As in `_reach`: a branch that could still reach an
                    # unmeasurable edge is never pruned.
                    if not onward_unmeasurable[depth_left].get(successor, False) and (
                        abs(new_product) * onward_bounds[depth_left].get(successor, 1.0)
                        < threshold
                    ):
                        continue
                stack.append(
                    (successor, trail + [successor], new_product, step_signed, unknown)
                )
        paths.sort(key=lambda entry: -abs(entry["product"]))
        return paths, truncated

    def attribute(
        self,
        target,
        source,
        *,
        threshold: float = 1e-6,
        max_depth: int = 32,
        max_paths: int = 512,
    ) -> PathAttribution:
        """Decompose ``source``'s influence on ``target`` path by path.

        Each path's contribution is the chain-rule product of its edge weights;
        the signed sum over paths is the end-to-end sensitivity, which is where
        cancellation between two routes shows up as a total far below the
        largest single path.

        Args:
            target: Destination node (id, port, or fragment).
            source: Origin node.
            threshold: Prune a path once ``|product|`` falls below this.
            max_depth: Maximum path length.
            max_paths: Stop after this many paths and mark the result
                truncated, rather than enumerating a combinatorial blow-up.
        """
        src = self.resolve(source)
        dst = self.resolve(target)
        paths, truncated = self._enumerate_paths(
            src, dst, threshold, max_depth, max_paths
        )
        all_signed = bool(paths) and all(
            entry["signed"] and not entry["unknown"] for entry in paths
        )
        return PathAttribution(
            target=dst,
            source=src,
            paths=paths,
            total=sum(entry["product"] for entry in paths) if all_signed else None,
            total_magnitude=sum(abs(entry["product"]) for entry in paths),
            truncated=truncated,
        )

    def dominant_paths(
        self,
        target,
        k: int = 5,
        *,
        source=None,
        threshold: float = 1e-6,
        max_depth: int = 32,
        max_paths: int = 512,
    ) -> List[Dict[str, Any]]:
        """The ``k`` strongest paths into ``target`` (optionally from ``source``).

        With no ``source``, every node with no in-edges inside the search — the
        model's genuine independent inputs and states — is used as an origin.
        """
        dst = self.resolve(target)
        if source is not None:
            return self.attribute(
                dst,
                source,
                threshold=threshold,
                max_depth=max_depth,
                max_paths=max_paths,
            ).paths[:k]

        # Reuse the slice rather than an unpruned reachability sweep: it applies
        # the same admissible bound, so the candidate origins are found without
        # an exponential search that could silently truncate.
        origin_slice = self.slice(dst, threshold, max_depth=max_depth)
        if origin_slice.truncated:
            warnings.warn(
                f"dominant_paths({target!r}): the slice used to find candidate "
                f"origins hit its search budget, so the ranking may be missing "
                f"stronger paths. Raise the threshold or lower max_depth.",
                UserWarning,
                stacklevel=2,
            )
        reachable = origin_slice.scores
        origins = [
            node
            for node in reachable
            if node != dst and self.graph.in_degree(node) == 0
        ]
        if not origins:
            # Every candidate origin is driven by something — a closed loop.
            # The states are then the model's independent variables.
            origins = [
                node
                for node in reachable
                if node != dst and self.graph.nodes[node]["kind"] == "state"
            ]
        collected: List[Dict[str, Any]] = []
        for origin in origins:
            paths, _ = self._enumerate_paths(
                origin, dst, threshold, max_depth, max_paths
            )
            collected.extend(paths)
        collected.sort(key=lambda entry: -abs(entry["product"]))
        return collected[:k]

    def dead_edges(self, threshold: float = 0.0) -> List[Dict[str, Any]]:
        """Structural edges that transmit no influence at this operating point.

        A wire the model declares and the mathematics ignores: a gain of zero, a
        saturated nonlinearity, a term that cancels. This is the quantitative
        form of a dead-store warning — the connection is real, the influence is
        not. Edges with no local gradient are excluded (unknown is not dead), and
        so are the state self-loops, whose zero A block is the *definition* of a
        plain integrator rather than a defect.
        """
        found = []
        for src, dst, data in self.graph.edges(data=True):
            if not data["local_gradient"] or src == dst:
                continue
            if abs(data["magnitude"]) <= threshold:
                found.append(
                    {
                        "src": src,
                        "dst": dst,
                        "kind": data["kind"],
                        "magnitude": abs(data["magnitude"]),
                    }
                )
        found.sort(key=lambda entry: (entry["magnitude"], entry["src"]))
        return found

    def bottlenecks(
        self,
        target,
        *,
        threshold: float = 0.01,
        max_depth: int = 32,
    ) -> List[str]:
        """Nodes every influential path to ``target`` must pass through.

        Computed on the slice at ``threshold``: a node is a bottleneck when
        deleting it disconnects at least one slice origin from ``target``. These
        are the signals worth instrumenting, and the single points of failure in
        a redundancy argument.

        Returns a bare list, so it has nowhere to report that the underlying
        slice was truncated — a truncated slice is missing paths, and a missing
        path is exactly what turns a non-bottleneck into an apparent one. That
        case warns instead; take the slice yourself and check
        :attr:`InfluenceSlice.truncated` if you need to handle it.
        """
        model_slice = self.slice(target, threshold, max_depth=max_depth)
        if model_slice.truncated:
            warnings.warn(
                f"bottlenecks({target!r}): the underlying slice hit its search "
                f"budget, so paths are missing and a node can look like a "
                f"single point of failure when it is not. Raise the threshold "
                f"or lower max_depth.",
                UserWarning,
                stacklevel=2,
            )
        subgraph = model_slice.subgraph
        subgraph.remove_edges_from(nx.selfloop_edges(subgraph))
        node = model_slice.target
        if node not in subgraph:
            return []
        origins = [
            other
            for other in subgraph
            if other != node and subgraph.in_degree(other) == 0
        ]
        if not origins:
            return []
        bottleneck = []
        for candidate in subgraph:
            if candidate == node or candidate in origins:
                continue
            trimmed = subgraph.copy()
            trimmed.remove_node(candidate)
            for origin in origins:
                if origin in trimmed and not nx.has_path(trimmed, origin, node):
                    bottleneck.append(candidate)
                    break
        return sorted(bottleneck)

    # -- reporting ---------------------------------------------------------

    def summary(self) -> str:
        """Human-readable overview: size, conventions, and honesty labels."""
        unknown = [
            (src, dst)
            for src, dst, data in self.graph.edges(data=True)
            if not data["local_gradient"]
        ]
        hybrid = sorted(
            {d["block"] for _, d in self.graph.nodes(data=True) if d["hybrid"]}
        )
        lines = [
            f"InfluenceGraph for {self.system.name}",
            f"  {self.n_blocks} blocks, {self.graph.number_of_nodes()} nodes, "
            f"{self.graph.number_of_edges()} edges",
            f"  at={self.at}  normalize={self.normalize}  tau={self.tau:g} s"
            f"  scale_floor={self.scale_floor:g}",
        ]
        if self.at == "trajectory":
            lines.append(
                f"  {len(self.times)} snapshots over "
                f"[{self.times[0]:g}, {self.times[-1]:g}] s, reduce={self.reduce}"
            )
        if unknown:
            lines.append(f"  {len(unknown)} edges with no local gradient:")
            for src, dst in unknown[:8]:
                note = self.graph.edges[src, dst]["note"]
                lines.append(f"    {src} -> {dst}  ({note})")
            if len(unknown) > 8:
                lines.append(f"    ... {len(unknown) - 8} more")
        if hybrid:
            lines.append(
                f"  {len(hybrid)} hybrid blocks (weights valid for the current "
                f"mode only): {', '.join(hybrid[:6])}"
            )
        floored = self.nodes_at_scale_floor()
        if floored:
            lines.append(
                f"  {len(floored)} signals sit at the scale floor (value ~0 at this "
                f"operating point), so elasticities through them are inflated by "
                f"the floor rather than measured: {', '.join(sorted(floored)[:5])}"
                + (f", ... (+{len(floored) - 5})" if len(floored) > 5 else "")
            )
            lines.append(
                "    Analyze at a settled operating point, raise scale_floor, or "
                "use normalize='none' if this matters."
            )
        dead = self.dead_edges()
        if dead:
            lines.append(f"  {len(dead)} dead edges (structural, zero influence)")
        return "\n".join(lines)

    def nodes_at_scale_floor(self) -> List[str]:
        """Signals whose normalizer came from ``scale_floor``, not from a value.

        A relative weight divides by the signal's magnitude, so a signal that is
        (near) zero at the operating point — an error signal at equilibrium, an
        integrator state at ``t=0`` — produces an elasticity governed by
        ``scale_floor`` rather than by the model. Those weights are not wrong so
        much as meaningless, and they are large, so they dominate any ranking.
        """
        if self.normalize != "relative":
            return []
        return sorted(
            node
            for node, data in self.graph.nodes(data=True)
            if "scale" in data
            and bool(np.any(np.asarray(data["scale"]) <= self.scale_floor))
        )

attribute(target, source, *, threshold=1e-06, max_depth=32, max_paths=512)

Decompose source's influence on target path by path.

Each path's contribution is the chain-rule product of its edge weights; the signed sum over paths is the end-to-end sensitivity, which is where cancellation between two routes shows up as a total far below the largest single path.

Parameters:

Name	Type	Description	Default
target		Destination node (id, port, or fragment).	required
source		Origin node.	required
threshold	float	Prune a path once |product| falls below this.	1e-06
max_depth	int	Maximum path length.	32
max_paths	int	Stop after this many paths and mark the result truncated, rather than enumerating a combinatorial blow-up.	512

Source code in jaxonomy/analysis/influence.py

def attribute(
    self,
    target,
    source,
    *,
    threshold: float = 1e-6,
    max_depth: int = 32,
    max_paths: int = 512,
) -> PathAttribution:
    """Decompose ``source``'s influence on ``target`` path by path.

    Each path's contribution is the chain-rule product of its edge weights;
    the signed sum over paths is the end-to-end sensitivity, which is where
    cancellation between two routes shows up as a total far below the
    largest single path.

    Args:
        target: Destination node (id, port, or fragment).
        source: Origin node.
        threshold: Prune a path once ``|product|`` falls below this.
        max_depth: Maximum path length.
        max_paths: Stop after this many paths and mark the result
            truncated, rather than enumerating a combinatorial blow-up.
    """
    src = self.resolve(source)
    dst = self.resolve(target)
    paths, truncated = self._enumerate_paths(
        src, dst, threshold, max_depth, max_paths
    )
    all_signed = bool(paths) and all(
        entry["signed"] and not entry["unknown"] for entry in paths
    )
    return PathAttribution(
        target=dst,
        source=src,
        paths=paths,
        total=sum(entry["product"] for entry in paths) if all_signed else None,
        total_magnitude=sum(abs(entry["product"]) for entry in paths),
        truncated=truncated,
    )

bottlenecks(target, *, threshold=0.01, max_depth=32)

Nodes every influential path to target must pass through.

Computed on the slice at threshold: a node is a bottleneck when deleting it disconnects at least one slice origin from target. These are the signals worth instrumenting, and the single points of failure in a redundancy argument.

Returns a bare list, so it has nowhere to report that the underlying slice was truncated — a truncated slice is missing paths, and a missing path is exactly what turns a non-bottleneck into an apparent one. That case warns instead; take the slice yourself and check :attr:InfluenceSlice.truncated if you need to handle it.

Source code in jaxonomy/analysis/influence.py

def bottlenecks(
    self,
    target,
    *,
    threshold: float = 0.01,
    max_depth: int = 32,
) -> List[str]:
    """Nodes every influential path to ``target`` must pass through.

    Computed on the slice at ``threshold``: a node is a bottleneck when
    deleting it disconnects at least one slice origin from ``target``. These
    are the signals worth instrumenting, and the single points of failure in
    a redundancy argument.

    Returns a bare list, so it has nowhere to report that the underlying
    slice was truncated — a truncated slice is missing paths, and a missing
    path is exactly what turns a non-bottleneck into an apparent one. That
    case warns instead; take the slice yourself and check
    :attr:`InfluenceSlice.truncated` if you need to handle it.
    """
    model_slice = self.slice(target, threshold, max_depth=max_depth)
    if model_slice.truncated:
        warnings.warn(
            f"bottlenecks({target!r}): the underlying slice hit its search "
            f"budget, so paths are missing and a node can look like a "
            f"single point of failure when it is not. Raise the threshold "
            f"or lower max_depth.",
            UserWarning,
            stacklevel=2,
        )
    subgraph = model_slice.subgraph
    subgraph.remove_edges_from(nx.selfloop_edges(subgraph))
    node = model_slice.target
    if node not in subgraph:
        return []
    origins = [
        other
        for other in subgraph
        if other != node and subgraph.in_degree(other) == 0
    ]
    if not origins:
        return []
    bottleneck = []
    for candidate in subgraph:
        if candidate == node or candidate in origins:
            continue
        trimmed = subgraph.copy()
        trimmed.remove_node(candidate)
        for origin in origins:
            if origin in trimmed and not nx.has_path(trimmed, origin, node):
                bottleneck.append(candidate)
                break
    return sorted(bottleneck)

dead_edges(threshold=0.0)

Structural edges that transmit no influence at this operating point.

A wire the model declares and the mathematics ignores: a gain of zero, a saturated nonlinearity, a term that cancels. This is the quantitative form of a dead-store warning — the connection is real, the influence is not. Edges with no local gradient are excluded (unknown is not dead), and so are the state self-loops, whose zero A block is the definition of a plain integrator rather than a defect.

Source code in jaxonomy/analysis/influence.py

def dead_edges(self, threshold: float = 0.0) -> List[Dict[str, Any]]:
    """Structural edges that transmit no influence at this operating point.

    A wire the model declares and the mathematics ignores: a gain of zero, a
    saturated nonlinearity, a term that cancels. This is the quantitative
    form of a dead-store warning — the connection is real, the influence is
    not. Edges with no local gradient are excluded (unknown is not dead), and
    so are the state self-loops, whose zero A block is the *definition* of a
    plain integrator rather than a defect.
    """
    found = []
    for src, dst, data in self.graph.edges(data=True):
        if not data["local_gradient"] or src == dst:
            continue
        if abs(data["magnitude"]) <= threshold:
            found.append(
                {
                    "src": src,
                    "dst": dst,
                    "kind": data["kind"],
                    "magnitude": abs(data["magnitude"]),
                }
            )
    found.sort(key=lambda entry: (entry["magnitude"], entry["src"]))
    return found

dominant_paths(target, k=5, *, source=None, threshold=1e-06, max_depth=32, max_paths=512)

The k strongest paths into target (optionally from source).

With no source, every node with no in-edges inside the search — the model's genuine independent inputs and states — is used as an origin.

Source code in jaxonomy/analysis/influence.py

def dominant_paths(
    self,
    target,
    k: int = 5,
    *,
    source=None,
    threshold: float = 1e-6,
    max_depth: int = 32,
    max_paths: int = 512,
) -> List[Dict[str, Any]]:
    """The ``k`` strongest paths into ``target`` (optionally from ``source``).

    With no ``source``, every node with no in-edges inside the search — the
    model's genuine independent inputs and states — is used as an origin.
    """
    dst = self.resolve(target)
    if source is not None:
        return self.attribute(
            dst,
            source,
            threshold=threshold,
            max_depth=max_depth,
            max_paths=max_paths,
        ).paths[:k]

    # Reuse the slice rather than an unpruned reachability sweep: it applies
    # the same admissible bound, so the candidate origins are found without
    # an exponential search that could silently truncate.
    origin_slice = self.slice(dst, threshold, max_depth=max_depth)
    if origin_slice.truncated:
        warnings.warn(
            f"dominant_paths({target!r}): the slice used to find candidate "
            f"origins hit its search budget, so the ranking may be missing "
            f"stronger paths. Raise the threshold or lower max_depth.",
            UserWarning,
            stacklevel=2,
        )
    reachable = origin_slice.scores
    origins = [
        node
        for node in reachable
        if node != dst and self.graph.in_degree(node) == 0
    ]
    if not origins:
        # Every candidate origin is driven by something — a closed loop.
        # The states are then the model's independent variables.
        origins = [
            node
            for node in reachable
            if node != dst and self.graph.nodes[node]["kind"] == "state"
        ]
    collected: List[Dict[str, Any]] = []
    for origin in origins:
        paths, _ = self._enumerate_paths(
            origin, dst, threshold, max_depth, max_paths
        )
        collected.extend(paths)
    collected.sort(key=lambda entry: -abs(entry["product"]))
    return collected[:k]

nodes_at_scale_floor()

Signals whose normalizer came from scale_floor, not from a value.

A relative weight divides by the signal's magnitude, so a signal that is (near) zero at the operating point — an error signal at equilibrium, an integrator state at t=0 — produces an elasticity governed by scale_floor rather than by the model. Those weights are not wrong so much as meaningless, and they are large, so they dominate any ranking.

Source code in jaxonomy/analysis/influence.py

def nodes_at_scale_floor(self) -> List[str]:
    """Signals whose normalizer came from ``scale_floor``, not from a value.

    A relative weight divides by the signal's magnitude, so a signal that is
    (near) zero at the operating point — an error signal at equilibrium, an
    integrator state at ``t=0`` — produces an elasticity governed by
    ``scale_floor`` rather than by the model. Those weights are not wrong so
    much as meaningless, and they are large, so they dominate any ranking.
    """
    if self.normalize != "relative":
        return []
    return sorted(
        node
        for node, data in self.graph.nodes(data=True)
        if "scale" in data
        and bool(np.any(np.asarray(data["scale"]) <= self.scale_floor))
    )

relative_threshold(target, fraction=0.01, *, direction='backward', max_depth=32, floor=1e-12)

A threshold set at fraction of the strongest influence on target.

An absolute threshold only reads as a percentage when tau is comparable to the time constants on the paths involved (see the module docstring). Scaling to the strongest score makes "keep what carries at least 1% of what the dominant contributor carries" mean the same thing at any tau.

Parameters:

Name	Type	Description	Default
target		Node id, port object, or name fragment.	required
fraction	float	Fraction of the strongest score to keep.	0.01
direction	str	As in :meth:slice.	'backward'
max_depth	int	As in :meth:slice.	32
floor	float	Threshold the reference sweep runs at, and the value returned when nothing upstream carries influence. It is passed to the search rather than left at zero so the sweep stays pruned; a model whose strongest contributor falls below it would yield floor itself.	1e-12

Returns:

Type	Description
float	A threshold to pass to :meth:slice / :meth:bottlenecks.

Source code in jaxonomy/analysis/influence.py

def relative_threshold(
    self,
    target,
    fraction: float = 0.01,
    *,
    direction: str = "backward",
    max_depth: int = 32,
    floor: float = 1e-12,
) -> float:
    """A threshold set at ``fraction`` of the strongest influence on ``target``.

    An absolute threshold only reads as a percentage when ``tau`` is
    comparable to the time constants on the paths involved (see the module
    docstring). Scaling to the strongest score makes "keep what carries at
    least 1% of what the dominant contributor carries" mean the same thing at
    any ``tau``.

    Args:
        target: Node id, port object, or name fragment.
        fraction: Fraction of the strongest score to keep.
        direction: As in :meth:`slice`.
        max_depth: As in :meth:`slice`.
        floor: Threshold the reference sweep runs at, and the value returned
            when nothing upstream carries influence. It is passed to the
            search rather than left at zero so the sweep stays pruned; a
            model whose strongest contributor falls below it would yield
            ``floor`` itself.

    Returns:
        A threshold to pass to :meth:`slice` / :meth:`bottlenecks`.
    """
    node = self.resolve(target)
    best, _, _, _ = self._reach(node, max_depth, direction, threshold=floor)
    others = [value for other, value in best.items() if other != node]
    if not others or max(others) <= 0.0:
        return floor
    return max(fraction * max(others), floor)

resolve(spec)

Turn a port object, locator, or name fragment into a node id.

Accepts an exact node id, an InputPort / OutputPort, a (system, port_index) locator, or any unambiguous suffix of a node id ("integ:out:out_0", "integ", "out:y").

Source code in jaxonomy/analysis/influence.py

def resolve(self, spec) -> str:
    """Turn a port object, locator, or name fragment into a node id.

    Accepts an exact node id, an ``InputPort`` / ``OutputPort``, a
    ``(system, port_index)`` locator, or any unambiguous suffix of a node
    id (``"integ:out:out_0"``, ``"integ"``, ``"out:y"``).
    """
    if isinstance(spec, (InputPort, OutputPort)):
        return self._require(port_node_id(spec))
    if isinstance(spec, tuple) and len(spec) == 2:
        system, index = spec
        ports = system.output_ports if hasattr(system, "output_ports") else []
        if index < len(ports):
            return self._require(port_node_id(ports[index]))
    if not isinstance(spec, str):
        raise TypeError(
            f"Cannot resolve {spec!r} to an influence-graph node; pass a node "
            f"id, a port object, or a (system, port_index) locator."
        )
    if spec in self.graph:
        return spec
    matches = [n for n in self.graph if n.endswith(spec) or spec in n]
    if len(matches) == 1:
        return matches[0]
    if not matches:
        raise KeyError(
            f"No influence-graph node matches {spec!r}. "
            f"Nodes are named '<block path>:in|out:<port>' or "
            f"'<block path>:xc|xd'; {self.graph.number_of_nodes()} exist, "
            f"e.g. {sorted(self.graph)[:4]}."
        )
    raise KeyError(
        f"{spec!r} is ambiguous — it matches {len(matches)} nodes: "
        f"{sorted(matches)[:6]}. Use a longer fragment or the full node id."
    )

slice(target, threshold=0.01, *, direction='backward', max_depth=32)

Quantitative model slice: what influences target by ≥ threshold.

The boolean answer — everything structurally upstream — is :meth:structural_slice; this one keeps only what lies on a path carrying at least threshold of the influence, in the relative-sensitivity sense described in the module docstring.

0.01 reads as "1%" only when tau is comparable to the time constants on the paths involved — the threshold is absolute, and a path across k integrators carries a factor of tau**k, so the same cutoff means different things at different tau. When the strongest contributor scores 95, threshold=0.01 retains everything down to ~0.01% of it, not 1%. Use :meth:relative_threshold to get the cutoff that means a fraction of the dominant contributor.

Two kinds of node are kept, and the distinction is load-bearing. A node is influential when its own best path to target clears the threshold. It is a connector when it merely lies on some influential node's best route: a relative weight is an elasticity, so a signal can pass through a junction that nearly cancels it and be amplified back afterwards, leaving a mid-route node with a small score of its own. Keeping only the influential ones would punch holes in the result — naming a block as influential while the route from it to the target ran through blocks that had been dropped, leaving :attr:InfluenceSlice.subgraph disconnected and :meth:bottlenecks meaningless. Connectors are read off the routes the search actually found, so nothing is added that no real path uses.

scores reports every retained node's own best product to the target, which is the number to rank by.

Parameters:

Name	Type	Description	Default
target		Node id, port object, or name fragment (see :meth:resolve).	required
threshold	float	Minimum |path product| for a path to be retained.	0.01
direction	str	"backward" (default, what influences the target) or "forward" (what the target influences).	'backward'
max_depth	int	Hard bound on path length, and the only hard bound — a partial product is not a bound on the whole path's (see :meth:_reach), so nothing may be pruned on the running value.	32

Returns:

Name	Type	Description
An	InfluenceSlice	class:InfluenceSlice.

Source code in jaxonomy/analysis/influence.py

def slice(
    self,
    target,
    threshold: float = 0.01,
    *,
    direction: str = "backward",
    max_depth: int = 32,
) -> InfluenceSlice:
    """Quantitative model slice: what influences ``target`` by ≥ ``threshold``.

    The boolean answer — everything structurally upstream — is
    :meth:`structural_slice`; this one keeps only what lies on a path
    carrying at least ``threshold`` of the influence, in the
    relative-sensitivity sense described in the module docstring.

    ``0.01`` reads as "1%" only when ``tau`` is comparable to the time
    constants on the paths involved — the threshold is absolute, and a
    path across *k* integrators carries a factor of ``tau**k``, so the
    same cutoff means different things at different ``tau``. When the
    strongest contributor scores 95, ``threshold=0.01`` retains everything
    down to ~0.01% of it, not 1%. Use :meth:`relative_threshold` to get
    the cutoff that means a fraction *of the dominant contributor*.

    Two kinds of node are kept, and the distinction is load-bearing. A node
    is **influential** when its own best path to ``target`` clears the
    threshold. It is a **connector** when it merely lies on some influential
    node's best route: a relative weight is an elasticity, so a signal can
    pass through a junction that nearly cancels it and be amplified back
    afterwards, leaving a mid-route node with a small score of its own.
    Keeping only the influential ones would punch holes in the result —
    naming a block as influential while the route from it to the target ran
    through blocks that had been dropped, leaving
    :attr:`InfluenceSlice.subgraph` disconnected and :meth:`bottlenecks`
    meaningless. Connectors are read off the routes the search actually
    found, so nothing is added that no real path uses.

    ``scores`` reports every retained node's own best product to the target,
    which is the number to rank by.

    Args:
        target: Node id, port object, or name fragment (see :meth:`resolve`).
        threshold: Minimum ``|path product|`` for a path to be retained.
        direction: ``"backward"`` (default, what influences the target) or
            ``"forward"`` (what the target influences).
        max_depth: Hard bound on path length, and the only hard bound — a
            partial product is not a bound on the whole path's (see
            :meth:`_reach`), so nothing may be pruned on the running value.

    Returns:
        An :class:`InfluenceSlice`.
    """
    if direction not in ("backward", "forward"):
        raise ValueError(
            f"direction must be 'backward' or 'forward', got {direction!r}"
        )
    node = self.resolve(target)
    best, unknown, routes, truncated = self._reach(
        node, max_depth, direction, threshold=threshold
    )

    influential = {
        other
        for other, score in best.items()
        if other == node or score >= threshold or unknown.get(other, False)
    }
    # Pull in whatever each influential node's own best route passes
    # through, so the result is a connected sub-model rather than a set of
    # names with no way to get between them.
    retained_nodes = set(influential)
    for other in influential:
        retained_nodes.update(routes.get(other, (other,)))

    # Every node on a retained route gets a score. Dropping the ones the
    # numeric search never scored would leave `blocks` naming fewer blocks
    # than `edges` actually connects — the retained set and the subgraph
    # would disagree, and the route's interior blocks would vanish from the
    # answer even though they are the only way the influence travels.
    retained = {
        other: best[other] if other in best else 0.0 for other in retained_nodes
    }
    retained.setdefault(node, best[node])

    kept = set()
    for other in influential:
        route = routes.get(other, ())
        for first, second in zip(route, route[1:]):
            kept.add(
                (first, second) if direction == "backward" else (second, first)
            )
    # Beyond the best routes, keep any edge between retained nodes that
    # itself carries a qualifying path — a parallel branch of comparable
    # strength belongs in the slice even though some other route was best.
    for other in retained:
        neighbours = (
            self.graph.successors(other)
            if direction == "backward"
            else self.graph.predecessors(other)
        )
        for downstream in neighbours:
            if downstream == other or downstream not in retained:
                continue
            edge = (
                (other, downstream)
                if direction == "backward"
                else (downstream, other)
            )
            magnitude, known = self._step_weight(*edge)
            if not known or magnitude * retained[downstream] >= threshold:
                kept.add(edge)

    blocks = sorted({self.graph.nodes[n]["block"] for n in retained})
    return InfluenceSlice(
        target=node,
        threshold=threshold,
        direction=direction,
        scores=dict(retained),
        edges=sorted(kept),
        blocks=blocks,
        unknown_nodes=sorted(
            other for other in retained if unknown.get(other, False)
        ),
        graph=self,
        truncated=truncated,
    )

structural_slice(target, *, direction='backward')

Boolean slice: every block structurally connected to target.

The over-approximation :meth:slice improves on, computed from the model's declared connectivity rather than from the weighted graph — so it stays a genuine bound even where a Jacobian could not be taken. Provided so the two can be compared directly on a real model.

Source code in jaxonomy/analysis/influence.py

def structural_slice(self, target, *, direction: str = "backward") -> List[str]:
    """Boolean slice: every block structurally connected to ``target``.

    The over-approximation :meth:`slice` improves on, computed from the
    model's declared connectivity rather than from the weighted graph — so
    it stays a genuine bound even where a Jacobian could not be taken.
    Provided so the two can be compared directly on a real model.
    """
    node = self.resolve(target)
    structure = self.structure if self.structure is not None else self.graph
    if node not in structure:
        return [self.graph.nodes[node]["block"]]
    reached = (
        nx.ancestors(structure, node)
        if direction == "backward"
        else nx.descendants(structure, node)
    )
    return sorted({structure.nodes[n]["block"] for n in reached | {node}})

summary()

Human-readable overview: size, conventions, and honesty labels.

Source code in jaxonomy/analysis/influence.py

def summary(self) -> str:
    """Human-readable overview: size, conventions, and honesty labels."""
    unknown = [
        (src, dst)
        for src, dst, data in self.graph.edges(data=True)
        if not data["local_gradient"]
    ]
    hybrid = sorted(
        {d["block"] for _, d in self.graph.nodes(data=True) if d["hybrid"]}
    )
    lines = [
        f"InfluenceGraph for {self.system.name}",
        f"  {self.n_blocks} blocks, {self.graph.number_of_nodes()} nodes, "
        f"{self.graph.number_of_edges()} edges",
        f"  at={self.at}  normalize={self.normalize}  tau={self.tau:g} s"
        f"  scale_floor={self.scale_floor:g}",
    ]
    if self.at == "trajectory":
        lines.append(
            f"  {len(self.times)} snapshots over "
            f"[{self.times[0]:g}, {self.times[-1]:g}] s, reduce={self.reduce}"
        )
    if unknown:
        lines.append(f"  {len(unknown)} edges with no local gradient:")
        for src, dst in unknown[:8]:
            note = self.graph.edges[src, dst]["note"]
            lines.append(f"    {src} -> {dst}  ({note})")
        if len(unknown) > 8:
            lines.append(f"    ... {len(unknown) - 8} more")
    if hybrid:
        lines.append(
            f"  {len(hybrid)} hybrid blocks (weights valid for the current "
            f"mode only): {', '.join(hybrid[:6])}"
        )
    floored = self.nodes_at_scale_floor()
    if floored:
        lines.append(
            f"  {len(floored)} signals sit at the scale floor (value ~0 at this "
            f"operating point), so elasticities through them are inflated by "
            f"the floor rather than measured: {', '.join(sorted(floored)[:5])}"
            + (f", ... (+{len(floored) - 5})" if len(floored) > 5 else "")
        )
        lines.append(
            "    Analyze at a settled operating point, raise scale_floor, or "
            "use normalize='none' if this matters."
        )
    dead = self.dead_edges()
    if dead:
        lines.append(f"  {len(dead)} dead edges (structural, zero influence)")
    return "\n".join(lines)

InfluenceSlice dataclass

A quantitative model slice: what actually reaches a target.

Attributes:

Name	Type	Description
target	str	Node id the slice was taken to (or from).
threshold	float	Influence cutoff a path had to clear to be included.
direction	str	"backward" (what influences the target) or "forward" (what the target influences).
scores	Dict[str, float]	{node_id: best |path product| between node and target}. For a node reached only across an edge with no local gradient the value is a bound, not a measurement — see unknown_nodes. Nodes kept only to connect an influential node to the target (see :meth:InfluenceGraph.slice) appear here with their own, possibly small, score.
edges	List[Tuple[str, str]]	(src, dst) pairs retained.
blocks	List[str]	Block name paths touched — the block-level slice.
unknown_nodes	List[str]	Nodes that some retained path reaches across an edge with no local gradient. Their score accounts only for the measurable routes, so it is not the whole story — treat it as a partial reading rather than a measurement.
truncated	bool	True if the path search hit its expansion budget, in which case the scores are lower bounds and the slice may be missing contributors.
graph	'InfluenceGraph'	The originating :class:InfluenceGraph.

Source code in jaxonomy/analysis/influence.py

@dataclass
class InfluenceSlice:
    """A quantitative model slice: what actually reaches a target.

    Attributes:
        target: Node id the slice was taken to (or from).
        threshold: Influence cutoff a path had to clear to be included.
        direction: ``"backward"`` (what influences the target) or
            ``"forward"`` (what the target influences).
        scores: ``{node_id: best |path product| between node and target}``. For
            a node reached only across an edge with no local gradient the value
            is a bound, not a measurement — see ``unknown_nodes``. Nodes kept
            only to connect an influential node to the target (see
            :meth:`InfluenceGraph.slice`) appear here with their own, possibly
            small, score.
        edges: ``(src, dst)`` pairs retained.
        blocks: Block name paths touched — the block-level slice.
        unknown_nodes: Nodes that some retained path reaches across an edge
            with no local gradient. Their score accounts only for the
            measurable routes, so it is not the whole story — treat it as a
            partial reading rather than a measurement.
        truncated: True if the path search hit its expansion budget, in which
            case the scores are lower bounds and the slice may be missing
            contributors.
        graph: The originating :class:`InfluenceGraph`.
    """

    target: str
    threshold: float
    direction: str
    scores: Dict[str, float]
    edges: List[Tuple[str, str]]
    blocks: List[str]
    unknown_nodes: List[str]
    graph: "InfluenceGraph"
    truncated: bool = False

    @property
    def unknown_paths(self) -> bool:
        """True if any retained path crosses an edge with no local gradient."""
        return bool(self.unknown_nodes)

    @property
    def block_scores(self) -> Dict[str, float]:
        """``{block name path: score}``, ranked, for the block-level answer.

        :attr:`scores` is keyed by *signal* (one node per input port, output
        port, and state group), which is the right granularity for tracing a
        route but the wrong one for "which block matters most". This reduces a
        block's nodes to one number by taking the **maximum**, so a block's
        score is that of its most influential signal.

        Max is the reducer because a block's input and output nodes lie on the
        *same* path — summing them would count one route twice, and a block's
        influence is not the sum of its ports' influences. The trade-off is
        that a block reached by several genuinely independent routes reads as
        its strongest one, not their total; use :meth:`attribute` when the
        split between routes is the question.

        The dict is ordered by descending score. Blocks holding a node in
        :attr:`unknown_nodes` are present with a score covering only their
        measurable routes — check that list before reading a rank as complete.
        """
        by_block: Dict[str, float] = {}
        for node, score in self.scores.items():
            block = self.graph.graph.nodes[node]["block"]
            by_block[block] = max(by_block.get(block, 0.0), abs(score))
        return dict(sorted(by_block.items(), key=lambda kv: (-kv[1], kv[0])))

    def __repr__(self) -> str:
        flags = ""
        if self.unknown_paths:
            flags += ", unknown paths"
        if self.truncated:
            flags += ", TRUNCATED"
        return (
            f"InfluenceSlice({self.target}, threshold={self.threshold:g}, "
            f"{len(self.blocks)} blocks, {len(self.scores)} nodes{flags})"
        )

    @property
    def subgraph(self) -> nx.DiGraph:
        """The retained portion of the influence graph.

        Built from the retained nodes *and* edges rather than as an edge-induced
        view, so a node with no retained edge — the target of a slice that keeps
        nothing else — is still present.
        """
        view = nx.DiGraph()
        for node in self.scores:
            view.add_node(node, **self.graph.graph.nodes[node])
        for src, dst in self.edges:
            view.add_edge(src, dst, **self.graph.graph.edges[src, dst])
        return view

    def report(self, by: str = "node") -> str:
        """Human-readable ranking.

        Args:
            by: ``"node"`` (default) ranks individual signals; ``"block"``
                ranks blocks via :attr:`block_scores`.
        """
        if by not in ("node", "block"):
            raise ValueError(f"by must be 'node' or 'block', got {by!r}")
        lines = [
            f"Influence slice ({self.direction}) for {self.target}",
            f"  threshold={self.threshold:g}  normalize={self.graph.normalize}  "
            f"tau={self.graph.tau:g}  at={self.graph.at}",
            f"  {len(self.blocks)} of {self.graph.n_blocks} blocks retained",
        ]
        unknown = set(self.unknown_nodes)
        if by == "block":
            unknown_blocks = {self.graph.graph.nodes[n]["block"] for n in unknown}
            # Drop the target *node*, not its whole block: another signal on the
            # same block (its input port feeding the target state, say) is a
            # real contributor and belongs in the ranking.
            ranked_blocks: Dict[str, float] = {}
            for node, score in self.scores.items():
                if node == self.target:
                    continue
                block = self.graph.graph.nodes[node]["block"]
                ranked_blocks[block] = max(ranked_blocks.get(block, 0.0), abs(score))
            for block, score in sorted(
                ranked_blocks.items(), key=lambda kv: (-kv[1], kv[0])
            ):
                flag = "  (bound)" if block in unknown_blocks else ""
                lines.append(f"    {score:>10.4g}  {block}{flag}")
        else:
            ranked = sorted(self.scores.items(), key=lambda kv: -abs(kv[1]))
            for node, score in ranked:
                if node == self.target:
                    continue
                flag = "  (bound)" if node in unknown else ""
                lines.append(f"    {score:>10.4g}  {node}{flag}")
        if unknown:
            lines.append(
                f"  NOTE: {len(unknown)} nodes are reached only across an edge with "
                f"no local gradient; they are kept unconditionally and their "
                f"scores are upper bounds."
            )
        if self.truncated:
            lines.append(
                "  NOTE: the path search hit its expansion budget — these scores "
                "are lower bounds and contributors may be missing. Lower "
                "max_depth or raise the threshold."
            )
        return "\n".join(lines)

block_scores property

{block name path: score}, ranked, for the block-level answer.

:attr:scores is keyed by signal (one node per input port, output port, and state group), which is the right granularity for tracing a route but the wrong one for "which block matters most". This reduces a block's nodes to one number by taking the maximum, so a block's score is that of its most influential signal.

Max is the reducer because a block's input and output nodes lie on the same path — summing them would count one route twice, and a block's influence is not the sum of its ports' influences. The trade-off is that a block reached by several genuinely independent routes reads as its strongest one, not their total; use :meth:attribute when the split between routes is the question.

The dict is ordered by descending score. Blocks holding a node in :attr:unknown_nodes are present with a score covering only their measurable routes — check that list before reading a rank as complete.

subgraph property

The retained portion of the influence graph.

Built from the retained nodes and edges rather than as an edge-induced view, so a node with no retained edge — the target of a slice that keeps nothing else — is still present.

unknown_paths property

True if any retained path crosses an edge with no local gradient.

report(by='node')

Human-readable ranking.

Parameters:

Name	Type	Description	Default
by	str	"node" (default) ranks individual signals; "block" ranks blocks via :attr:block_scores.	'node'

Source code in jaxonomy/analysis/influence.py

def report(self, by: str = "node") -> str:
    """Human-readable ranking.

    Args:
        by: ``"node"`` (default) ranks individual signals; ``"block"``
            ranks blocks via :attr:`block_scores`.
    """
    if by not in ("node", "block"):
        raise ValueError(f"by must be 'node' or 'block', got {by!r}")
    lines = [
        f"Influence slice ({self.direction}) for {self.target}",
        f"  threshold={self.threshold:g}  normalize={self.graph.normalize}  "
        f"tau={self.graph.tau:g}  at={self.graph.at}",
        f"  {len(self.blocks)} of {self.graph.n_blocks} blocks retained",
    ]
    unknown = set(self.unknown_nodes)
    if by == "block":
        unknown_blocks = {self.graph.graph.nodes[n]["block"] for n in unknown}
        # Drop the target *node*, not its whole block: another signal on the
        # same block (its input port feeding the target state, say) is a
        # real contributor and belongs in the ranking.
        ranked_blocks: Dict[str, float] = {}
        for node, score in self.scores.items():
            if node == self.target:
                continue
            block = self.graph.graph.nodes[node]["block"]
            ranked_blocks[block] = max(ranked_blocks.get(block, 0.0), abs(score))
        for block, score in sorted(
            ranked_blocks.items(), key=lambda kv: (-kv[1], kv[0])
        ):
            flag = "  (bound)" if block in unknown_blocks else ""
            lines.append(f"    {score:>10.4g}  {block}{flag}")
    else:
        ranked = sorted(self.scores.items(), key=lambda kv: -abs(kv[1]))
        for node, score in ranked:
            if node == self.target:
                continue
            flag = "  (bound)" if node in unknown else ""
            lines.append(f"    {score:>10.4g}  {node}{flag}")
    if unknown:
        lines.append(
            f"  NOTE: {len(unknown)} nodes are reached only across an edge with "
            f"no local gradient; they are kept unconditionally and their "
            f"scores are upper bounds."
        )
    if self.truncated:
        lines.append(
            "  NOTE: the path search hit its expansion budget — these scores "
            "are lower bounds and contributors may be missing. Lower "
            "max_depth or raise the threshold."
        )
    return "\n".join(lines)

LeafJacobians dataclass

Local Jacobian blocks for one leaf at one operating point.

Attributes:

Name	Type	Description
leaf	Any	The LeafSystem these Jacobians describe.
u0	List[Any]	Operating-point value of each input port, in port order.
y0	List[Any]	Operating-point value of each output port, in port order.
x0	Dict[str, Any]	Operating-point value of each state kind present, keyed by "xc" / "xd".
d	Dict[Tuple[int, int], ndarray]	{(out_i, in_j): ndarray(m_i, n_j)} — direct feedthrough.
c	Dict[Tuple[str, int], ndarray]	{(kind, out_i): ndarray(m_i, n_x)} — state → output.
b	Dict[Tuple[str, int], ndarray]	{(kind, in_j): ndarray(n_x, n_j)} — input → state rate/update.
a	Dict[Tuple[str, str], ndarray]	{(src_kind, dst_kind): ndarray(n_dst, n_src)} — state → state rate/update, including the cross terms (an ODE reading discrete state, a periodic update reading continuous state).
notes	Dict[str, str]	{subject: reason} for every quantity that could not be differentiated, e.g. {"out:mode": "non-inexact dtype int32"}. Callers turn these into local_gradient=None edge labels rather than silently reporting a zero.

Source code in jaxonomy/analysis/block_jacobians.py

@dataclass
class LeafJacobians:
    """Local Jacobian blocks for one leaf at one operating point.

    Attributes:
        leaf: The ``LeafSystem`` these Jacobians describe.
        u0: Operating-point value of each input port, in port order.
        y0: Operating-point value of each output port, in port order.
        x0: Operating-point value of each state kind present, keyed by
            ``"xc"`` / ``"xd"``.
        d: ``{(out_i, in_j): ndarray(m_i, n_j)}`` — direct feedthrough.
        c: ``{(kind, out_i): ndarray(m_i, n_x)}`` — state → output.
        b: ``{(kind, in_j): ndarray(n_x, n_j)}`` — input → state rate/update.
        a: ``{(src_kind, dst_kind): ndarray(n_dst, n_src)}`` — state → state
            rate/update, including the cross terms (an ODE reading discrete
            state, a periodic update reading continuous state).
        notes: ``{subject: reason}`` for every quantity that could *not* be
            differentiated, e.g. ``{"out:mode": "non-inexact dtype int32"}``.
            Callers turn these into ``local_gradient=None`` edge labels rather
            than silently reporting a zero.
    """

    leaf: Any
    u0: List[Any]
    y0: List[Any]
    x0: Dict[str, Any]
    d: Dict[Tuple[int, int], np.ndarray] = field(default_factory=dict)
    c: Dict[Tuple[str, int], np.ndarray] = field(default_factory=dict)
    b: Dict[Tuple[str, int], np.ndarray] = field(default_factory=dict)
    a: Dict[Tuple[str, str], np.ndarray] = field(default_factory=dict)
    notes: Dict[str, str] = field(default_factory=dict)

    def __repr__(self) -> str:
        return (
            f"LeafJacobians({self.leaf.name}, n_in={len(self.u0)}, "
            f"n_out={len(self.y0)}, states={sorted(self.x0)}, "
            f"notes={len(self.notes)})"
        )

PathAttribution dataclass

Chain-rule decomposition of one source's influence on one target.

Attributes:

Name	Type	Description
target	str	Destination node id.
source	str	Origin node id.
paths	List[Dict[str, Any]]	One entry per path, ranked by |product|, each a dict with nodes, product, signed (False when a matrix block on the path made the sign meaningless) and unknown (True when an edge on the path has no local gradient).
total	Optional[float]	Signed sum of path products when every path is signed, else None — a sum of magnitudes would hide cancellation.
total_magnitude	float	Sum of |product| over paths, always available.
truncated	bool	True if enumeration hit max_paths or max_depth.

Source code in jaxonomy/analysis/influence.py

@dataclass
class PathAttribution:
    """Chain-rule decomposition of one source's influence on one target.

    Attributes:
        target: Destination node id.
        source: Origin node id.
        paths: One entry per path, ranked by ``|product|``, each a dict with
            ``nodes``, ``product``, ``signed`` (False when a matrix block on
            the path made the sign meaningless) and ``unknown`` (True when an
            edge on the path has no local gradient).
        total: Signed sum of path products when every path is signed, else
            ``None`` — a sum of magnitudes would hide cancellation.
        total_magnitude: Sum of ``|product|`` over paths, always available.
        truncated: True if enumeration hit ``max_paths`` or ``max_depth``.
    """

    target: str
    source: str
    paths: List[Dict[str, Any]]
    total: Optional[float]
    total_magnitude: float
    truncated: bool

    def __repr__(self) -> str:
        total = "n/a" if self.total is None else f"{self.total:.4g}"
        return (
            f"PathAttribution({self.source} -> {self.target}, "
            f"{len(self.paths)} paths, total={total})"
        )

    def report(self, max_paths: int = 10) -> str:
        lines = [f"Attribution {self.source} -> {self.target}"]
        if self.total is not None:
            lines.append(f"  total (signed sum over paths): {self.total:.6g}")
        lines.append(f"  total magnitude: {self.total_magnitude:.6g}")
        for entry in self.paths[:max_paths]:
            flags = []
            if not entry["signed"]:
                flags.append("magnitude-only")
            if entry["unknown"]:
                flags.append("no local gradient on path")
            suffix = f"  [{', '.join(flags)}]" if flags else ""
            lines.append(f"  {entry['product']:+.6g}{suffix}")
            lines.append(f"      {' -> '.join(entry['nodes'])}")
        if len(self.paths) > max_paths:
            lines.append(f"  ... {len(self.paths) - max_paths} more paths")
        if self.truncated:
            lines.append("  NOTE: enumeration truncated; totals are partial.")
        return "\n".join(lines)

format_influence_subgraph(graph, focus, edges, types=None, rates=None, dropped_for_budget=0)

Render a node/edge selection as compact, citable text.

The footer distinguishes the two reasons a block can be absent, because to a reader they mean opposite things. A block left out because its influence fell below the threshold is known to be negligible — that is an answer. A block left out because the budget ran out is simply unknown, and treating it as negligible would be a fabrication. Without the footer both look identical: missing.

Source code in jaxonomy/analysis/influence_context.py

def format_influence_subgraph(
    graph: InfluenceGraph,
    focus: List[str],
    edges: List[Tuple[str, str]],
    types: Optional[Dict[str, str]] = None,
    rates: Optional[Dict[str, str]] = None,
    dropped_for_budget: int = 0,
) -> str:
    """Render a node/edge selection as compact, citable text.

    The footer distinguishes the two reasons a block can be absent, because to
    a reader they mean opposite things. A block left out because its influence
    fell below the threshold is *known to be negligible* — that is an answer. A
    block left out because the budget ran out is simply *unknown*, and treating
    it as negligible would be a fabrication. Without the footer both look
    identical: missing.
    """
    types = types if types is not None else _block_types(graph)
    rates = rates if rates is not None else _block_rates(graph)
    nodes = list(dict.fromkeys(list(focus) + [n for edge in edges for n in edge]))

    lines = [
        f"# influence subgraph of model '{graph.system.name}'",
        f"# focus: {', '.join(focus)}",
        f"# weights: {'relative sensitivity (dimensionless)' if graph.normalize == 'relative' else 'raw partial derivative'}"
        f"; state-rate edges scaled by tau={graph.tau:g} s"
        f"; evaluated at={graph.at}"
        + (f" (reduce={graph.reduce})" if graph.at == "trajectory" else ""),
        "# signals: id | value | units | sample_time | block_type",
    ]
    for node in nodes:
        data = graph.graph.nodes[node]
        block = data["block"]
        marker = " *" if node in focus else ""
        lines.append(
            f"{node}{marker} | {_value_text(data['value'])} | "
            f"{_units_text(data['units'])} | {rates.get(block, '?')} | "
            f"{types.get(block, '?')}"
            + (" | hybrid" if data["hybrid"] else "")
        )
    lines.append("# edges: src -> dst | weight | kind")
    for src, dst in edges:
        data = graph.graph.edges[src, dst]
        line = f"{src} -> {dst} | {_weight_text(data)} | {data['kind']}"
        if data["note"]:
            line += f" | {data['note']}"
        lines.append(line)

    shown = {graph.graph.nodes[n]["block"] for n in nodes}
    negligible = graph.n_blocks - len(shown)
    lines.append(
        f"# coverage: {len(shown)}/{graph.n_blocks} blocks; {negligible} omitted "
        f"as below-cutoff (negligible, NOT unknown)"
    )
    if dropped_for_budget:
        lines.append(
            f"# WARNING: {dropped_for_budget} edges also dropped for budget — "
            f"anything reachable only via those is UNKNOWN, not negligible"
        )
    return "\n".join(lines)

influence_graph(system, context=None, *, at='operating_point', results=None, times=None, n_snapshots=5, tau=1.0, normalize='relative', scale_floor=1e-06, probe=None, reduce='max', simulator_options=None)

Build the sensitivity-weighted influence graph of a model.

Parameters:

Name	Type	Description	Default
system		A Diagram or a single LeafSystem.	required
context		Root context fixing the operating point. Defaults to system.create_context().	None
at	str	"operating_point" weights every edge once, at context. "trajectory" weights at several snapshots and stores per-edge profiles — the honest answer when a nonlinearity means one number per edge cannot be right everywhere (a block saturated at the operating point has a zero local gradient there and a large one elsewhere).	'operating_point'
results		A SimulationResults supplying the snapshot times for at="trajectory". The states are re-derived by advancing context, because recorded signals do not pin down every stateful leaf — which costs one simulate call per snapshot. Budget for that on a large model: simulate's fixed setup cost scales with block count and dominates the integration itself (a 1 µs span costs the same as a 4 s one), so n_snapshots=6 on a 2500-block model is minutes rather than seconds. Building at a single operating point is linear in block count and stays in seconds at that size.	None
times	Optional[Sequence[float]]	Explicit snapshot times, used instead of results.	None
n_snapshots	int	How many times to take from results.time.	5
tau	float	Seconds of integration represented by a continuous-state-rate edge; only affects edges into ẋc. Set it from the fastest state on the paths you care about — every integrator on a path contributes a factor of tau, so a value taken from the slow dynamics of a stiff model inflates multi-integrator path products (see the module docstring).	1.0
normalize	str	"relative" (default, dimensionless elasticities) or "none" (raw partial derivatives in model units).	'relative'
scale_floor	float	Floor on a signal's operating-point magnitude when normalizing, so a signal that happens to sit at zero does not produce an infinite elasticity. Nodes at the floor are visible via their value attribute.	1e-06
probe	Optional[float]	When set to a relative step size (0.05 = 5% of each signal's magnitude), every edge whose exact derivative is zero is re-checked with a central-difference secant, and the secant is used instead when it is non-zero. This is the cross-check for the one thing an exact local derivative gets wrong: a quantizer between steps, a saturation at its rail, or a dead zone inside the zone is locally flat while still transmitting information, and would otherwise be reported dead. Costs two extra block evaluations per signal component; None (default) skips it.	None
reduce	str	How a trajectory profile collapses to the scalar weight used by queries: "max" (default, conservative — never hides an influence that appears at some point), "mean", or "final".	'max'
simulator_options		SimulatorOptions for the trajectory-mode re-integration.	None

Returns:

Name	Type	Description
An	InfluenceGraph	class:InfluenceGraph.

Example

import jaxonomy from jaxonomy.library import Constant, Gain, Integrator from jaxonomy.analysis import influence_graph builder = jaxonomy.DiagramBuilder() source = builder.add(Constant(1.0, name="src")) gain = builder.add(Gain(3.0, name="gain")) plant = builder.add(Integrator(1.0, name="plant")) builder.connect(source.output_ports[0], gain.input_ports[0]) builder.connect(gain.output_ports[0], plant.input_ports[0]) diagram = builder.build(name="root") graph = influence_graph(diagram) graph.slice("plant:xc", threshold=0.01).blocks ['gain', 'plant', 'src']

Source code in jaxonomy/analysis/influence.py

def influence_graph(
    system,
    context=None,
    *,
    at: str = "operating_point",
    results=None,
    times: Optional[Sequence[float]] = None,
    n_snapshots: int = 5,
    tau: float = 1.0,
    normalize: str = "relative",
    scale_floor: float = 1e-6,
    probe: Optional[float] = None,
    reduce: str = "max",
    simulator_options=None,
) -> InfluenceGraph:
    """Build the sensitivity-weighted influence graph of a model.

    Args:
        system: A ``Diagram`` or a single ``LeafSystem``.
        context: Root context fixing the operating point. Defaults to
            ``system.create_context()``.
        at: ``"operating_point"`` weights every edge once, at ``context``.
            ``"trajectory"`` weights at several snapshots and stores per-edge
            profiles — the honest answer when a nonlinearity means one number
            per edge cannot be right everywhere (a block saturated at the
            operating point has a zero local gradient there and a large one
            elsewhere).
        results: A ``SimulationResults`` supplying the snapshot times for
            ``at="trajectory"``. The states are re-derived by advancing
            ``context``, because recorded signals do not pin down every
            stateful leaf — which costs one ``simulate`` call per snapshot.
            Budget for that on a large model: ``simulate``'s fixed setup cost
            scales with block count and dominates the integration itself (a
            1 µs span costs the same as a 4 s one), so ``n_snapshots=6`` on a
            2500-block model is minutes rather than seconds. Building at a
            single operating point is linear in block count and stays in
            seconds at that size.
        times: Explicit snapshot times, used instead of ``results``.
        n_snapshots: How many times to take from ``results.time``.
        tau: Seconds of integration represented by a continuous-state-rate
            edge; only affects edges into ``ẋc``. Set it from the *fastest*
            state on the paths you care about — every integrator on a path
            contributes a factor of ``tau``, so a value taken from the slow
            dynamics of a stiff model inflates multi-integrator path products
            (see the module docstring).
        normalize: ``"relative"`` (default, dimensionless elasticities) or
            ``"none"`` (raw partial derivatives in model units).
        scale_floor: Floor on a signal's operating-point magnitude when
            normalizing, so a signal that happens to sit at zero does not
            produce an infinite elasticity. Nodes at the floor are visible via
            their ``value`` attribute.
        probe: When set to a relative step size (``0.05`` = 5% of each signal's
            magnitude), every edge whose exact derivative is zero is re-checked
            with a central-difference secant, and the secant is used instead when
            it is non-zero. This is the cross-check for the one thing an exact
            local derivative gets wrong: a quantizer between steps, a saturation
            at its rail, or a dead zone inside the zone is *locally* flat while
            still transmitting information, and would otherwise be reported dead.
            Costs two extra block evaluations per signal component; ``None``
            (default) skips it.
        reduce: How a trajectory profile collapses to the scalar weight used by
            queries: ``"max"`` (default, conservative — never hides an
            influence that appears at some point), ``"mean"``, or ``"final"``.
        simulator_options: ``SimulatorOptions`` for the trajectory-mode
            re-integration.

    Returns:
        An :class:`InfluenceGraph`.

    Example:
        >>> import jaxonomy
        >>> from jaxonomy.library import Constant, Gain, Integrator
        >>> from jaxonomy.analysis import influence_graph
        >>> builder = jaxonomy.DiagramBuilder()
        >>> source = builder.add(Constant(1.0, name="src"))
        >>> gain = builder.add(Gain(3.0, name="gain"))
        >>> plant = builder.add(Integrator(1.0, name="plant"))
        >>> builder.connect(source.output_ports[0], gain.input_ports[0])
        >>> builder.connect(gain.output_ports[0], plant.input_ports[0])
        >>> diagram = builder.build(name="root")
        >>> graph = influence_graph(diagram)
        >>> graph.slice("plant:xc", threshold=0.01).blocks
        ['gain', 'plant', 'src']
    """
    if at not in ("operating_point", "trajectory"):
        raise ValueError(
            f"at must be 'operating_point' or 'trajectory', got {at!r}"
        )
    if normalize not in ("relative", "none"):
        raise ValueError(
            f"normalize must be 'relative' or 'none', got {normalize!r}"
        )
    if reduce not in ("max", "mean", "final"):
        raise ValueError(f"reduce must be 'max', 'mean' or 'final', got {reduce!r}")
    if tau <= 0:
        raise ValueError(f"tau must be positive (seconds), got {tau}")

    if context is None:
        context = system.create_context()

    if at == "operating_point":
        jacobians, secants, notes = _local_data(system, context, probe, scale_floor)
        graph = _build_at(
            system,
            jacobians,
            tau=tau,
            normalize=normalize,
            scale_floor=scale_floor,
            secants=secants,
            probe=probe,
        )
        return InfluenceGraph(
            system=system,
            graph=graph,
            tau=tau,
            normalize=normalize,
            scale_floor=scale_floor,
            at=at,
            times=None,
            reduce=reduce,
            block_notes=notes,
            structure=_structural_graph(system),
        )

    snapshot_times = _snapshot_times(context, results, times, n_snapshots)
    contexts = _trajectory_contexts(system, context, snapshot_times, simulator_options)

    per_snapshot = []
    merged_notes: Dict[str, Dict[str, str]] = {}
    for snapshot_context in contexts:
        jacobians, secants, notes = _local_data(
            system, snapshot_context, probe, scale_floor
        )
        per_snapshot.append((jacobians, secants))
        for block, block_notes in notes.items():
            merged_notes.setdefault(block, {}).update(block_notes)

    scales = _trajectory_scales(
        system, [jacobians for jacobians, _ in per_snapshot], normalize, scale_floor
    )
    graphs = [
        _build_at(
            system,
            jacobians,
            tau=tau,
            normalize=normalize,
            scale_floor=scale_floor,
            scales=scales,
            secants=secants,
            probe=probe,
        )
        for jacobians, secants in per_snapshot
    ]

    combined = _merge_profiles(graphs, reduce)
    return InfluenceGraph(
        system=system,
        graph=combined,
        tau=tau,
        normalize=normalize,
        scale_floor=scale_floor,
        at=at,
        times=snapshot_times,
        reduce=reduce,
        block_notes=merged_notes,
        structure=_structural_graph(system),
    )

influence_subgraph(graph, focus, *, budget_tokens=1500, hops=4, threshold=0.0, direction='both')

A bounded, budgeted, citable neighbourhood of focus.

Parameters:

Name	Type	Description	Default
graph	InfluenceGraph	An :class:~jaxonomy.analysis.influence.InfluenceGraph.	required
focus		One or more focus points — node ids, port objects, name fragments, or a block name (which expands to all of that block's signals).	required
budget_tokens	int	Approximate ceiling on the rendered text, at :data:CHARS_PER_TOKEN characters per token. Edges are dropped weakest-first to fit; the result reports what was dropped.	1500
hops	int	How many graph edges out from the focus to expand. Nodes are signals, so crossing one block costs two hops (wire in, block Jacobian out) — the default of 4 reaches roughly two blocks.	4
threshold	float	Minimum edge |weight| to include at all.	0.0
direction	str	"both" (default), "backward" (what influences the focus) or "forward" (what it influences).	'both'

Returns:

Type	Description
Dict[str, Any]	A dict with text (the rendered context), nodes, edges,
Dict[str, Any]	blocks, focus, estimated_tokens, dropped_edges and
Dict[str, Any]	conventions. The dict is JSON-serializable.

Source code in jaxonomy/analysis/influence_context.py

def influence_subgraph(
    graph: InfluenceGraph,
    focus,
    *,
    budget_tokens: int = 1500,
    hops: int = 4,
    threshold: float = 0.0,
    direction: str = "both",
) -> Dict[str, Any]:
    """A bounded, budgeted, citable neighbourhood of ``focus``.

    Args:
        graph: An :class:`~jaxonomy.analysis.influence.InfluenceGraph`.
        focus: One or more focus points — node ids, port objects, name
            fragments, or a block name (which expands to all of that block's
            signals).
        budget_tokens: Approximate ceiling on the rendered text, at
            :data:`CHARS_PER_TOKEN` characters per token. Edges are dropped
            weakest-first to fit; the result reports what was dropped.
        hops: How many graph edges out from the focus to expand. Nodes are
            signals, so crossing one block costs two hops (wire in, block
            Jacobian out) — the default of 4 reaches roughly two blocks.
        threshold: Minimum edge ``|weight|`` to include at all.
        direction: ``"both"`` (default), ``"backward"`` (what influences the
            focus) or ``"forward"`` (what it influences).

    Returns:
        A dict with ``text`` (the rendered context), ``nodes``, ``edges``,
        ``blocks``, ``focus``, ``estimated_tokens``, ``dropped_edges`` and
        ``conventions``. The dict is JSON-serializable.
    """
    if direction not in ("both", "backward", "forward"):
        raise ValueError(
            f"direction must be 'both', 'backward' or 'forward', got {direction!r}"
        )
    if budget_tokens <= 0:
        raise ValueError(f"budget_tokens must be positive, got {budget_tokens}")

    focus_nodes = _resolve_focus(graph, focus)
    ranked = _expand(graph, focus_nodes, hops, threshold, direction)

    types = _block_types(graph)
    rates = _block_rates(graph)

    # Grow the edge set strongest-first until the rendered text would exceed
    # the budget.  Rendering is cheap relative to building the graph, so the
    # budget is enforced on the real text rather than on an estimate of it.
    kept: List[Tuple[str, str]] = []
    dropped: List[Tuple[str, str]] = []
    text = format_influence_subgraph(graph, focus_nodes, [], types, rates)
    for index, (_hop, _negative, edge) in enumerate(ranked):
        candidate = kept + [edge]
        remaining = len(ranked) - index - 1
        rendered = format_influence_subgraph(
            graph, focus_nodes, candidate, types, rates, dropped_for_budget=remaining
        )
        if len(rendered) > budget_tokens * CHARS_PER_TOKEN and kept:
            # Edges arrive strongest-first, so once one does not fit, none of
            # the weaker ones would have been a better use of the budget.
            dropped = [entry[2] for entry in ranked[index:]]
            break
        kept, text = candidate, rendered
    # Re-render with the final count so the warning matches what was dropped.
    text = format_influence_subgraph(
        graph, focus_nodes, kept, types, rates, dropped_for_budget=len(dropped)
    )

    nodes = list(dict.fromkeys(focus_nodes + [n for edge in kept for n in edge]))
    return {
        "model": graph.system.name,
        "focus": focus_nodes,
        "text": text,
        "estimated_tokens": -(-len(text) // CHARS_PER_TOKEN),
        "blocks": sorted({graph.graph.nodes[n]["block"] for n in nodes}),
        "nodes": [
            {
                "id": node,
                "kind": graph.graph.nodes[node]["kind"],
                "block": graph.graph.nodes[node]["block"],
                "block_type": types.get(graph.graph.nodes[node]["block"], "?"),
                "port": graph.graph.nodes[node]["port"],
                "size": graph.graph.nodes[node]["size"],
                "value": _value_text(graph.graph.nodes[node]["value"]),
                "units": _units_text(graph.graph.nodes[node]["units"]),
                "sample_time": rates.get(graph.graph.nodes[node]["block"], "?"),
                "hybrid": bool(graph.graph.nodes[node]["hybrid"]),
            }
            for node in nodes
        ],
        "edges": [
            {
                "src": src,
                "dst": dst,
                "kind": graph.graph.edges[src, dst]["kind"],
                "weight": _weight_text(graph.graph.edges[src, dst]),
                "local_gradient": bool(
                    graph.graph.edges[src, dst]["local_gradient"]
                ),
                "note": graph.graph.edges[src, dst]["note"],
            }
            for src, dst in kept
        ],
        "dropped_edges": [{"src": src, "dst": dst} for src, dst in dropped],
        "conventions": {
            "at": graph.at,
            "normalize": graph.normalize,
            "tau_seconds": graph.tau,
            "scale_floor": graph.scale_floor,
            "reduce": graph.reduce if graph.at == "trajectory" else None,
            "weight_meaning": (
                "relative (elasticity) sensitivity; a path's product is the "
                "relative end-to-end sensitivity"
                if graph.normalize == "relative"
                else "raw partial derivative in model units"
            ),
        },
    }

leaf_jacobians(leaf, root_context)

Compute every local Jacobian block of leaf at root_context.

Parameters:

Name	Type	Description	Default
leaf		A LeafSystem belonging to the system root_context was created from.	required
root_context		Root context supplying the operating point — time, parameters, this leaf's state, and (via upstream evaluation) the values arriving on its input ports.	required

Returns:

Name	Type	Description
A	LeafJacobians	class:LeafJacobians. Blocks that could not be differentiated are
	LeafJacobians	absent from d / c / b / a and explained in notes;
	LeafJacobians	this function does not raise on a non-differentiable block.

Source code in jaxonomy/analysis/block_jacobians.py

def leaf_jacobians(leaf, root_context) -> LeafJacobians:
    """Compute every local Jacobian block of ``leaf`` at ``root_context``.

    Args:
        leaf: A ``LeafSystem`` belonging to the system ``root_context`` was
            created from.
        root_context: Root context supplying the operating point — time,
            parameters, this leaf's state, and (via upstream evaluation) the
            values arriving on its input ports.

    Returns:
        A :class:`LeafJacobians`.  Blocks that could not be differentiated are
        absent from ``d`` / ``c`` / ``b`` / ``a`` and explained in ``notes``;
        this function does not raise on a non-differentiable block.
    """
    # An enabled port cache would hand back stored constants instead of
    # re-evaluating under the JVP, zeroing every gradient.  Inputs are fixed
    # below, so nothing upstream needs the cache.
    if root_context.port_cache:
        root_context = root_context.with_port_cache({})

    leaf_context = root_context[leaf.system_id]
    try:
        u0 = [port.eval(root_context) for port in leaf.input_ports]
    except Exception as exc:  # noqa: BLE001 - e.g. a dangling, unfixed input
        return LeafJacobians(
            leaf=leaf,
            u0=[],
            y0=[],
            x0={},
            notes={
                "block": (
                    f"input-port evaluation failed: {type(exc).__name__}: {exc}"
                )
            },
        )

    x0: Dict[str, Any] = {}
    if leaf_context.has_continuous_state and leaf.ode_callback is not None:
        x0["xc"] = leaf_context.continuous_state
    if leaf_context.has_discrete_state and _discrete_update_event(leaf) is not None:
        x0["xd"] = leaf_context.discrete_state

    jacs = LeafJacobians(leaf=leaf, u0=u0, y0=[], x0=x0)

    if leaf_context.has_discrete_state and "xd" not in x0:
        jacs.notes["state:xd"] = (
            "discrete state with no single periodic update callback"
        )

    def evaluate(us, xc, xd):
        return _eval_leaf(leaf, root_context, us, xc, xd)

    # Baseline evaluation, and the reference for output sizes / dtypes.
    try:
        y0, xcdot0, xdp0 = evaluate(u0, None, None)
    except Exception as exc:  # noqa: BLE001 - a block we cannot evaluate at all
        jacs.notes["block"] = f"evaluation failed: {type(exc).__name__}: {exc}"
        return jacs
    jacs.y0 = y0

    # Which outputs / states are differentiable at all.
    out_ok: Dict[int, int] = {}
    for i, y in enumerate(y0):
        reason = _differentiable(y)
        if reason is None:
            out_ok[i] = _size(y)
        else:
            jacs.notes[f"out:{leaf.output_ports[i].name}"] = reason

    in_ok: Dict[int, int] = {}
    for j, u in enumerate(u0):
        reason = _differentiable(u)
        if reason is None:
            in_ok[j] = _size(u)
        else:
            jacs.notes[f"in:{leaf.input_ports[j].name}"] = reason

    rates = {}  # state kind -> (baseline rate/update value, its size)
    if "xc" in x0 and xcdot0 is not None:
        rates["xc"] = (xcdot0, _size(xcdot0))
    if "xd" in x0 and xdp0 is not None:
        rates["xd"] = (xdp0, _size(xdp0))

    state_ok: Dict[str, int] = {}
    for kind, value in x0.items():
        reason = _differentiable(value)
        if reason is None and kind in rates:
            state_ok[kind] = _size(value)
        elif reason is not None:
            jacs.notes[f"state:{kind}"] = reason

    # --- ∂(·)/∂uⱼ : one forward-mode trace per input port -------------------
    for j, n_j in in_ok.items():
        u_flat, unravel = ravel_pytree(u0[j])

        def wrt_input(u_vec, j=j, unravel=unravel):
            us = list(u0)
            us[j] = unravel(u_vec)
            outputs, xcdot, xd_plus = evaluate(us, None, None)
            return (
                [_ravel(outputs[i]) for i in sorted(out_ok)],
                _ravel(xcdot) if "xc" in rates else None,
                _ravel(xd_plus) if "xd" in rates else None,
            )

        try:
            jac_out, jac_xcdot, jac_xdp = jax.jacfwd(wrt_input)(u_flat)
        except Exception as exc:  # noqa: BLE001
            jacs.notes[f"in:{leaf.input_ports[j].name}"] = (
                f"not differentiable: {type(exc).__name__}: {exc}"
            )
            continue

        for slot, i in enumerate(sorted(out_ok)):
            jacs.d[(i, j)] = _as_matrix(jac_out[slot], out_ok[i], n_j)
        for kind, jac in (("xc", jac_xcdot), ("xd", jac_xdp)):
            if jac is not None and kind in rates:
                jacs.b[(kind, j)] = _as_matrix(jac, rates[kind][1], n_j)

    # --- ∂(·)/∂x : one forward-mode trace per state kind -------------------
    for kind, n_x in state_ok.items():
        x_flat, unravel = ravel_pytree(x0[kind])

        def wrt_state(x_vec, kind=kind, unravel=unravel):
            xc = unravel(x_vec) if kind == "xc" else None
            xd = unravel(x_vec) if kind == "xd" else None
            outputs, xcdot, xd_plus = evaluate(u0, xc, xd)
            return (
                [_ravel(outputs[i]) for i in sorted(out_ok)],
                _ravel(xcdot) if "xc" in rates else None,
                _ravel(xd_plus) if "xd" in rates else None,
            )

        try:
            jac_out, jac_xcdot, jac_xdp = jax.jacfwd(wrt_state)(x_flat)
        except Exception as exc:  # noqa: BLE001
            jacs.notes[f"state:{kind}"] = (
                f"not differentiable: {type(exc).__name__}: {exc}"
            )
            continue

        for slot, i in enumerate(sorted(out_ok)):
            jacs.c[(kind, i)] = _as_matrix(jac_out[slot], out_ok[i], n_x)
        for dst, jac in (("xc", jac_xcdot), ("xd", jac_xdp)):
            if jac is not None and dst in rates:
                jacs.a[(kind, dst)] = _as_matrix(jac, rates[dst][1], n_x)

    return jacs