diff --git a/balm/_sd_algorithms/compute_attractor_seeds.py b/balm/_sd_algorithms/compute_attractor_seeds.py index 8b9564b..1bc9f9a 100644 --- a/balm/_sd_algorithms/compute_attractor_seeds.py +++ b/balm/_sd_algorithms/compute_attractor_seeds.py @@ -30,7 +30,7 @@ def compute_attractor_seeds( if balm.succession_diagram.DEBUG: print(f"[{node_id}] Start computing attractor seeds.") - node_space = sd.node_space(node_id) + node_space = sd.node_data(node_id)["space"] if len(node_space) == sd.network.variable_count(): # This node is a fixed-point. @@ -39,8 +39,8 @@ def compute_attractor_seeds( # Compute the list of child spaces if the node is expanded. Otherwise # "pretend" that there are no children. child_spaces = [] - if sd.node_is_expanded(node_id): - child_spaces = [sd.node_space(s) for s in sd.node_successors(node_id)] + if sd.node_data(node_id)["expanded"]: + child_spaces = [sd.node_data(s)["space"] for s in sd.node_successors(node_id)] # Fix everything in the NFVS to zero, as long as # it isn't already fixed by our `node_space`. diff --git a/balm/_sd_algorithms/expand_attractor_seeds.py b/balm/_sd_algorithms/expand_attractor_seeds.py index dfa4c79..1cd2f29 100644 --- a/balm/_sd_algorithms/expand_attractor_seeds.py +++ b/balm/_sd_algorithms/expand_attractor_seeds.py @@ -47,7 +47,7 @@ def expand_attractor_seeds(sd: SuccessionDiagram, size_limit: int | None = None) # Retrieve the stable motifs of children that are already expanded. expanded_children = [ - x for x in sd.node_successors(node) if sd.node_is_expanded(x) + x for x in sd.node_successors(node) if sd.node_data(x)["expanded"] ] expanded_motifs = [ sd.edge_stable_motif(node, child) for child in expanded_children @@ -61,7 +61,7 @@ def expand_attractor_seeds(sd: SuccessionDiagram, size_limit: int | None = None) # and continue to the next one. successors.pop() continue - if sd.node_is_expanded(successors[-1]): + if sd.node_data(successors[-1])["expanded"]: # The next node to explore is expanded (by some previous procedure) # but not "seen" in this search yet. We need to visit this node # regardless of other conditions @@ -69,7 +69,7 @@ def expand_attractor_seeds(sd: SuccessionDiagram, size_limit: int | None = None) # Now, we need to asses if the next successor has some candidate states which # are not covered by the already expanded children. - successor_space = sd.node_space(successors[-1]) + successor_space = sd.node_data(successors[-1])["space"] retained_set = make_retained_set( sd.symbolic, sd.node_nfvs(node), successor_space ) diff --git a/balm/_sd_algorithms/expand_minimal_spaces.py b/balm/_sd_algorithms/expand_minimal_spaces.py index e100dc1..761d589 100644 --- a/balm/_sd_algorithms/expand_minimal_spaces.py +++ b/balm/_sd_algorithms/expand_minimal_spaces.py @@ -34,7 +34,7 @@ def expand_minimal_spaces(sd: SuccessionDiagram, size_limit: int | None = None) successors = sorted(successors, reverse=True) # For determinism! # (reversed because we explore the list from the back) - node_space = sd.node_space(node) + node_space = sd.node_data(node)["space"] # Remove all immediate successors that are already visited or those who # do not cover any new minimal trap space. @@ -52,7 +52,7 @@ def expand_minimal_spaces(sd: SuccessionDiagram, size_limit: int | None = None) # of this node is already in the succession diagram. if len(successors) == 0: if sd.node_is_minimal(node): - minimal_traps.remove(sd.node_space(node)) + minimal_traps.remove(sd.node_data(node)["space"]) continue # At this point, we know that `s` is not visited and it contains diff --git a/balm/_sd_algorithms/expand_to_target.py b/balm/_sd_algorithms/expand_to_target.py index adb6ce9..1348172 100644 --- a/balm/_sd_algorithms/expand_to_target.py +++ b/balm/_sd_algorithms/expand_to_target.py @@ -25,7 +25,7 @@ def expand_to_target( while len(current_level) > 0: for node in current_level: - node_space = sd.node_space(node) + node_space = sd.node_data(node)["space"] if intersect(node_space, target) is None: # If `node_space` does not intersect with `target`, it is not relevant diff --git a/balm/control.py b/balm/control.py index 4603d1d..6a71358 100644 --- a/balm/control.py +++ b/balm/control.py @@ -58,7 +58,7 @@ def __init__( Example ------- >>> import balm - >>> sd = balm.SuccessionDiagram.from_bnet( + >>> sd = balm.SuccessionDiagram.from_rules( ... \"\"\" ... A, B & C ... B, A & C @@ -87,7 +87,7 @@ def __init__( self._control: list[ControlOverrides] = [] for c in control: cs = sorted(map(lambda x: sorted(x.items()), c)) - self._control.append(list(map(dict, cs))) + self._control.append(list(map(dict, cs))) # type: ignore self._strategy = strategy self._succession = succession @@ -198,7 +198,7 @@ def succession_control( ------- >>> import balm >>> from balm.control import succession_control - >>> sd = balm.SuccessionDiagram.from_bnet( + >>> sd = balm.SuccessionDiagram.from_rules( ... \"\"\" ... S, S ... A, S | B @@ -237,7 +237,7 @@ def succession_control( ------- >>> import balm >>> from balm.control import succession_control - >>> sd = balm.SuccessionDiagram.from_bnet( + >>> sd = balm.SuccessionDiagram.from_rules( ... \"\"\" ... S, S ... A, S | B @@ -335,7 +335,7 @@ def successions_to_target( ) for s in succession_diagram.node_ids(): - fixed_vars = succession_diagram.node_space(s) + fixed_vars = succession_diagram.node_data(s)["space"] if not is_subspace(fixed_vars, target): continue diff --git a/balm/drivers.py b/balm/drivers.py index 445626a..3cfeebc 100644 --- a/balm/drivers.py +++ b/balm/drivers.py @@ -40,7 +40,7 @@ def find_single_node_LDOIs( Example _______ >>> import balm - >>> sd = balm.SuccessionDiagram.from_bnet('A, A\\nB, A') + >>> sd = balm.SuccessionDiagram.from_rules('A, A\\nB, A') >>> ldois = balm.drivers.find_single_node_LDOIs(sd.network) >>> for k in sorted(ldois): ... print(f'{k} ==> {sorted(ldois[k].items())}') @@ -107,7 +107,7 @@ def find_single_drivers( Example ------- >>> import balm - >>> sd = balm.SuccessionDiagram.from_bnet('A, A\\nB, A') + >>> sd = balm.SuccessionDiagram.from_rules('A, A\\nB, A') >>> drivers = balm.drivers.find_single_drivers({'B': 0}, sd.network) >>> sorted(list(drivers)) [('A', 0), ('B', 0)] diff --git a/balm/interaction_graph_utils.py b/balm/interaction_graph_utils.py index da6dcc4..7b58cfd 100644 --- a/balm/interaction_graph_utils.py +++ b/balm/interaction_graph_utils.py @@ -96,7 +96,7 @@ def feedback_vertex_set( -------- >>> import balm >>> from balm.interaction_graph_utils import feedback_vertex_set - >>> sd = balm.SuccessionDiagram.from_bnet(\"\"\" + >>> sd = balm.SuccessionDiagram.from_rules(\"\"\" ... A, B ... B, A ... C, D diff --git a/balm/motif_avoidant.py b/balm/motif_avoidant.py index 13e2c69..4def880 100644 --- a/balm/motif_avoidant.py +++ b/balm/motif_avoidant.py @@ -1,5 +1,5 @@ """ - A module responsible for detecting motif-avoidant attractors within terminal restriction space. +A module responsible for detecting motif-avoidant attractors within terminal restriction space. """ from __future__ import annotations @@ -103,13 +103,16 @@ def detect_motif_avoidant_attractors( """ Compute a sub-list of `candidates` which correspond to motif-avoidant attractors. Other method inputs: - - `graph` and `petri_net` represent the model in which the property - should be checked. - - `terminal_restriction_space` is a symbolic set of states which contains - all motif avoidant attractors (i.e. if a candidate state can leave this - set, the candidate cannot be an attractor). - - `max_iterations` specifies how much time should be spent on the "simpler" - preprocessing before applying a more complete method. + + `graph` and `petri_net` represent the model in which the property + should be checked. + + `terminal_restriction_space` is a symbolic set of states which contains + all motif avoidant attractors (i.e. if a candidate state can leave this + set, the candidate cannot be an attractor). + + `max_iterations` specifies how much time should be spent on the "simpler" + preprocessing before applying a more complete method. """ if ensure_subspace is None: ensure_subspace = {} diff --git a/balm/succession_diagram.py b/balm/succession_diagram.py index 5d8fb80..d8875ca 100644 --- a/balm/succession_diagram.py +++ b/balm/succession_diagram.py @@ -1,6 +1,6 @@ from __future__ import annotations -from typing import TYPE_CHECKING, Any, cast +from typing import TYPE_CHECKING, Any, Literal, cast if TYPE_CHECKING: from typing import Iterator @@ -23,84 +23,63 @@ ) from balm.space_utils import percolate_space, space_unique_key from balm.trappist_core import trappist -from balm.types import BooleanSpace, SuccessionDiagramState +from balm.types import BooleanSpace, NodeData, SuccessionDiagramState -# Enables helpful "progress" messages. -DEBUG = False +DEBUG: bool = False +""" +If `True`, print debug and progress messages. +""" class SuccessionDiagram: """ - `SuccessionDiagram` (SD) is a directed acyclic graph representing the - structure of trap spaces induced by a particular `BooleanNetwork`. The root - of the graph is the whole state space (after percolation of constant - values). The leaf nodes are individual minimal trap spaces. Each path from - the root to a leaf represents a succession of gradually more restrictive - trap spaces. - - There are several additional "features" implemented by `SuccessionDiagram` - that allow us to implement more advanced (and efficient) algorithms: - - The succession diagram can be expanded lazily. Each node is initially - a *stub* node, - meaning none of its child nodes are known. A stub can be then *expanded* - into a full node by computing the *stable motifs* for the associated - trap space. - Each node can be annotated with *attractor seed states*, - i.e. states that are known to cover all network attractors within that - space. - - Overview of the `SuccessionDiagram` API: - - Introspection of the whole diagram: - * `SuccessionDiagram.root()`: the integer ID of the root node. - * `SuccessionDiagram.depth()`: The depth (hight?) of the current - succession diagram. - * `SuccessionDiagram.node_ids() / .stub_ids() / .expanded_ids()`: - Iterates over the - corresponding *node IDs* managed by this SD. * - `SuccessionDiagram.is_subgraph() / .is_isomorphic()`: Compares two - succession diagrams. - - Inspecting individual nodes: - * `SuccessionDiagram.node_depth(id)`: The depth (length of a maximal - path from the - root) of the given node. * `SuccessionDiagram.node_space(id)`: The - trap space associated with the given node. * - `SuccessionDiagram.node_is_minimal(id)`: Checks if the given node is - a minimal trap space: i.e. it is expanded and has no successor - nodes. * `SuccessionDiagram.node_is_expaned(id)`: Check if the given - node is expanded, i.e. its successor nodes are computed. * - `SuccessionDiagram.node_successors(id, compute=True/False)`: The - list of successor node IDs for the given node (can be computed if - not yet known, in which case the node becomes expanded). * - `SuccessionDiagram.node_attractor_seeds(id, compute=True/False)`: - The list of "attractor seed states" associated with the given node - (if these are not known yet, they can be computed). * - `SuccessionDiagram.edge_stable_motif(id, id)`: Obtain the stable - motif which enables the edge between the two nodes. - - There are several "expand procedures" that can explore a larger part - of the SD at once, - typically with a specific goal in mind: - * `SuccessionDiagram.expand_bfs() / .expand_dfs()`: Expands the - whole SD up to a certain - upper bound on size/depth. * - `SuccessionDiagram.expand_minimal_traps()`: Expands the SD such that - each minimal trap space is reachable by at least one path from the - root. - - *Other internal implementation notes:* - - *The node IDs are assumed to be a continuous range of integers. If this - breaks at any point, please make sure to check the correctness of the new - implementation thoroughly.* - - *There is no way to remove or otherwise "contract" the succession diagram: - Once a node is expanded, it stays expanded. At the time of writing, there - does not appear to be a need for deleting SD nodes and it greatly simplifies - reasoning about correctness.* + Succession diagram of a Boolean network. + + This encodes relationships between minimal trap spaces and can be used for + attractor detection and control. Note that the succession diagram is + expanded lazily, so it is not built until `build` or a similar method is + called. + + + Examples + -------- + >>> import balm + >>> sd = balm.SuccessionDiagram.from_rules(\""" + ... A, B + ... B, A & C + ... C, !A | B + ... \""") + >>> print(sd.summary()) # not built yet! + Succession Diagram with 1 nodes and depth 0. + State order: A, B, C + + Attractors in diagram: + + >>> sd.build() + >>> print(sd.summary()) # now it's built + Succession Diagram with 3 nodes and depth 1. + State order: A, B, C + + Attractors in diagram: + + minimal trap space 001 + ...................001 + + minimal trap space 111 + ...................111 + + """ + NFVS_NODE_THRESHOLD: int = 2_000 + """ + If the number of nodes in the (redunced) network is greater than this + threshold, we find a feedback vertex set (without consideration of cycle + sign) for reachability preprocessing. There is a trade-off between the speed + gains from a smaller node set to consider and the cost of determining which + FVS nodes only intersect negative cycles to find an NFVS subset. Typcially, + for smaller networks, the trade-off is in favor of computing a smaller NFVS. """ - NFVS_NODE_THRESHOLD = ( - 2_000 # if more than this many nodes, we only find FVS, not NFVS - ) __slots__ = ( "network", "symbolic", @@ -112,21 +91,44 @@ class SuccessionDiagram: def __init__(self, network: BooleanNetwork): # Original Boolean network. - self.network = cleanup_network(network) - self.symbolic = AsynchronousGraph(self.network) - # A Petri net representation of the original Boolean network. - self.petri_net = network_to_petrinet(network) + self.network: BooleanNetwork = cleanup_network(network) + """ + The Boolean network as an AEON network. + """ + + self.symbolic: AsynchronousGraph = AsynchronousGraph(self.network) + """ + The symbolic representation of the network. + """ + + self.petri_net: nx.DiGraph = network_to_petrinet(network) + """ + The Petri net representation of the network. + """ + if DEBUG: print( f"Generated global Petri net with {len(self.petri_net.nodes)} nodes and {len(self.petri_net.edges)} edges." ) + # Initially, we don't need the NFVS for anything. self.nfvs: list[str] | None = None - # A directed acyclic graph representing the succession diagram. - self.dag = nx.DiGraph() + """ + The negative feedback vertex set used for attractor detection. + """ + + self.dag: nx.DiGraph = nx.DiGraph() + """ + The directed acyclic graph (DAG) representation of the succession + diagram structure. + """ # A dictionary used for uniqueness checks on the nodes of the succession # diagram. See `SuccessionDiagram.ensure_node` for details. self.node_indices: dict[int, int] = {} + """ + A dictionary mapping subspace keys to their positions in the succession + diagram (see :func:`balm.space_utils.space_unique_key`). + """ # Create an un-expanded root node. self._ensure_node(None, {}) @@ -156,25 +158,35 @@ def __len__(self) -> int: return self.dag.number_of_nodes() @staticmethod - def from_aeon(model: str) -> SuccessionDiagram: - """ - Read a `BooleanNetwork` from the string contents of an `.aeon` model. - """ - return SuccessionDiagram(BooleanNetwork.from_aeon(model)) - - @staticmethod - def from_bnet(model: str) -> SuccessionDiagram: - """ - Read a `BooleanNetwork` from the string contents of a `.bnet` model. - """ - return SuccessionDiagram(BooleanNetwork.from_bnet(model)) - - @staticmethod - def from_sbml(model: str) -> SuccessionDiagram: - """ - Read a `BooleanNetwork` from the string contents of an `.sbml` model. - """ - return SuccessionDiagram(BooleanNetwork.from_sbml(model)) + def from_rules( + rules: str, + format: Literal["bnet", "aeon", "sbml"] = "bnet", + ) -> SuccessionDiagram: + """ + Generate a succession diagram from the given string. + + Parameters + ---------- + rules : str + The string representation of the network rules. + format : Literal['bnet', 'aeon', 'sbml'] + The format of the string. One of `"bnet"`, `"aeon"`, or `"sbml"`. + Defaults to `"bnet"`. + + Returns + ------- + SuccessionDiagram + The generated succession diagram. Initially unexpanded. + """ + + if format == "bnet": + return SuccessionDiagram(BooleanNetwork.from_bnet(rules)) + elif format == "aeon": + return SuccessionDiagram(BooleanNetwork.from_aeon(rules)) + elif format == "sbml": + return SuccessionDiagram(BooleanNetwork.from_sbml(rules)) + else: + raise ValueError(f"Unknown format: {format}") @staticmethod def from_file(path: str) -> SuccessionDiagram: @@ -184,12 +196,52 @@ def from_file(path: str) -> SuccessionDiagram: """ return SuccessionDiagram(BooleanNetwork.from_file(path)) - def expanded_attractor_seeds(self) -> list[list[BooleanSpace]]: - return [self.node_attractor_seeds(id) for id in self.expanded_ids()] + def expanded_attractor_seeds(self) -> dict[int, list[BooleanSpace]]: + """ + Attractor seeds for each expanded node. + + An attractor seed is a state belonging to an attractor (and thus a state + from which the entire attractor is reachable, by definition). + + If called before the `SuccessionDiagram` is fully built, this will not + be a complete accounting of attractor seed states. + + Returns + ------- + dict[int,list[BooleanSpace]] + Each key is the id of an expanded succession diagram node, whose + corresponding value is a list of attractor seeds for that node. Note + that one succession diagram node can have multiple attractors. Ordering + of the lists in the returned dictionary is not guaranteed. + + Example + ------- + >>> import balm + >>> sd = balm.SuccessionDiagram.from_rules(\""" + ... A, B + ... B, A & C + ... C, !A | B + ... \""") + >>> sd.build() + >>> eas = sd.expanded_attractor_seeds() + >>> for id, atts in sorted(eas.items()): + ... for x in atts: + ... print(f"{id}: {dict(sorted(x.items()))}") + 1: {'A': 0, 'B': 0, 'C': 1} + 2: {'A': 1, 'B': 1, 'C': 1} + """ + res: dict[int, list[BooleanSpace]] = {} + for id in self.expanded_ids(): + atts = self.node_attractor_seeds(id) + if not atts: # no attractors for this node + continue + res[id] = atts + + return res def summary(self) -> str: """ - Return a summary of the succession diagram. + Return a summary of the succession diagram as a string. """ var_ordering = sorted( [self.network.get_variable_name(v) for v in self.network.variables()] @@ -205,7 +257,10 @@ def summary(self) -> str: except KeyError: continue - space = self.node_space(node) + if not attrs: + continue + + space = self.node_data(node)["space"] if self.node_is_minimal(node): space_str_prefix = "minimal trap space " @@ -222,7 +277,8 @@ def summary(self) -> str: attr_str = "".join(str(v) for _, v in sorted(attr.items())) report_string += "." * len(space_str_prefix) + f"{attr_str}\n" report_string += "\n" - return report_string + # remove final extra newline and return + return report_string[:-1] def root(self) -> int: """ @@ -238,7 +294,7 @@ def depth(self) -> int: """ d = 0 for node in cast(set[int], self.dag.nodes()): - d = max(d, self.node_depth(int(node))) + d = max(d, self.node_data(int(node))["depth"]) return d def node_ids(self) -> Iterator[int]: @@ -253,7 +309,7 @@ def stub_ids(self) -> Iterator[int]: Iterator over all node IDs that are currently not expanded. """ for i in range(len(self)): - if not self.node_is_expanded(i): + if not self.node_data(i)["expanded"]: yield i def expanded_ids(self) -> Iterator[int]: @@ -261,7 +317,7 @@ def expanded_ids(self) -> Iterator[int]: Iterator over all node IDs that are currently expanded. """ for i in range(len(self)): - if self.node_is_expanded(i): + if self.node_data(i)["expanded"]: yield i def minimal_trap_spaces(self) -> list[int]: @@ -277,6 +333,17 @@ def find_node(self, node_space: BooleanSpace) -> int | None: """ Return the ID of the node matching the provided `node_space`, or `None` if no such node exists in this succession diagram. + + Parameters + ---------- + node_space: BooleanSpace + The space of the node to find. + + Returns + ------- + int | None + The ID of the node matching the provided `node_space`, or `None` + if no such node exists in this succession diagram. """ try: key = space_unique_key(node_space, self.network) # throws IndexError @@ -303,20 +370,31 @@ def is_subgraph(self, other: SuccessionDiagram) -> bool: WARNING: This does not take into account the stable motifs on individual edges. Just the subspaces associated with nodes and the presence of edges between nodes. + + Parameters + ---------- + other: SuccessionDiagram + The other succession diagram. + + Returns + ------- + bool + `True` if this succession diagram is a subgraph of the `other` + succession diagram. """ # Every stub node is reachable through an expanded node and # thus will be checked by the following code. for i in self.expanded_ids(): - other_i = other.find_node(self.node_space(i)) + other_i = other.find_node(self.node_data(i)["space"]) if other_i is None: return False my_successors = self.node_successors(i) other_successors = [] - if other.node_is_expanded(other_i): + if other.node_data(other_i)["expanded"]: other_successors = other.node_successors(other_i) for my_s in my_successors: - other_s = other.find_node(self.node_space(my_s)) + other_s = other.find_node(self.node_data(my_s)["space"]) if other_s not in other_successors: return False return True @@ -332,46 +410,73 @@ def is_isomorphic(self, other: SuccessionDiagram) -> bool: WARNING: This does not take into account the stable motifs on individual edges. Just the subspaces associated with nodes and the presence of edges between nodes. + + Parameters + ---------- + other: SuccessionDiagram + The other succession diagram. + + Returns + ------- + bool + `True` if the two succession diagrams are isomorphic. """ return self.is_subgraph(other) and other.is_subgraph(self) - def node_depth(self, node_id: int) -> int: - """ - Get the depth associated with the provided `node_id`. The depth is counted - as the longest path from the root node to the given node. + def node_data(self, node_id: int) -> NodeData: """ - return cast(int, self.dag.nodes[node_id]["depth"]) + Get the data associated with the provided `node_id`. - def node_space(self, node_id: int) -> BooleanSpace: - """ - Get the sub-space associated with the provided `node_id`. + Returns a `NodeData` object with the following attributes: - Note that this is the space *after* percolation. Hence it can hold that - `|node_space(child)| < |node_space(parent)| + |stable_motif(parent, child)|`. - """ - return cast(BooleanSpace, self.dag.nodes[node_id]["space"]) + - `depth`: The depth of the node. + - `attractors`: The attractors of the node. + - `petri_net`: The Petri net representation of the node. + - `space`: The sub-space of the node. + - `expanded`: Whether the node is expanded or not. - def node_is_expanded(self, node_id: int) -> bool: - """ - True if the successors of the given node are already computed. + See :class:`balm.types.NodeData` for more information. + + Parameters + ---------- + node_id: int + The ID of the node. + + Returns + ------- + NodeData + The data associated with the provided `node_id`. Note that at + runtime, this object is an untyped dictionary. """ - return cast(bool, self.dag.nodes[node_id]["expanded"]) + return cast(NodeData, self.dag.nodes[node_id]) def node_is_minimal(self, node_id: int) -> bool: """ - True if the node is expanded and it has no successors, i.e. it is a - minimal trap space. + True if the node represents a minimal trap space. + + Parameters + ---------- + node_id: int + The ID of the node. + + Returns + ------- + bool + `True` if the node is expanded and it has no successors, i.e. it is a + minimal trap space. """ is_leaf: bool = self.dag.out_degree(node_id) == 0 # type: ignore return is_leaf and self.dag.nodes[node_id]["expanded"] # type: ignore def node_successors(self, node_id: int, compute: bool = False) -> list[int]: """ - Return the successor nodes for the given `node_id`. If the node is - already expanded, known results are simply returned. If the node is not - expanded, but `compute` is set to `True`, then the node is expanded and - the newly computed results are returned. If the node is not expanded and - `compute` is set to `False`, the method raises a `KeyError` exception. + Return the successor nodes for the given `node_id`. + + If the node is already expanded, known results are simply returned. If + the node is not expanded, but `compute` is set to `True`, then the node + is expanded and the newly computed results are returned. If the node is + not expanded and `compute` is set to `False`, the method raises a + `KeyError` exception. The default behaviour intentionally does not compute successors to prevent "accidental complexity". @@ -383,6 +488,18 @@ def node_successors(self, node_id: int, compute: bool = False) -> list[int]: Also note that if the given `node_id` already has associated attractor data but is not expanded, this data will be deleted as it is no longer up to date. + + Parameters + ---------- + node_id: int + The ID of the node. + compute: bool + Whether to compute the successors if they are not already known. + + Returns + ------- + list[int] + The list of successor node ids. """ node = cast(dict[str, Any], self.dag.nodes[node_id]) if not node["expanded"] and not compute: @@ -396,8 +513,9 @@ def node_attractor_seeds( self, node_id: int, compute: bool = False ) -> list[BooleanSpace]: """ - Return the list of attractor seed states corresponding to the given - `node_id`. Similar to `node_successors`, the method either computes the + Return the list of attractor seed states for the given `node_id`. + + Similar to :meth:`node_successors`, the method either computes the data if unknown, or throws an exception, depending on the `compute` flag. @@ -405,10 +523,22 @@ def node_attractor_seeds( attractors are not guaranteed to be unique (i.e. you can "discover" the same attractor in multiple stub nodes, if the stub nodes intersect), and (b) this data is erased if the stub node is expanded later on. + + Parameters + ---------- + node_id: int + The ID of the node. + compute: bool + Whether to compute the attractor seeds if they are not already known. + + Returns + ------- + list[BooleanSpace] + The list of attractor seed states. """ - node = cast(dict[str, Any], self.dag.nodes[node_id]) + node = cast(NodeData, self.dag.nodes[node_id]) - attractors = cast(list[BooleanSpace] | None, node["attractors"]) + attractors = node["attractors"] if attractors is None and not compute: raise KeyError(f"Attractor data not computed for node {node_id}.") @@ -421,8 +551,19 @@ def node_attractor_seeds( def node_nfvs(self, node_id: int) -> list[str]: """ - Return the approximation of a minimum negative feedback vertex set that - is valid for the specified SD node. + Approximate minimum negative feedback vertex set on the subspace for the given node. + + See :func:`balm.interaction_graph_utils.feedback_vertex_set` for further details. + + Parameters + ---------- + node_id: int + The ID of the node. + + Returns + ------- + list[str] + The negative feedback vertex set, as a list of node names. """ # TODO: Right now, we are just returning the gloval FVS. But in the future, @@ -440,30 +581,33 @@ def node_nfvs(self, node_id: int) -> list[str]: return self.nfvs - def node_restricted_petri_net(self, node_id: int) -> nx.DiGraph | None: - """ - Return the pre-computed Petri net representation restricted to the subspace - of the specified SD node. - - This can return `None` if the requested node is already fully expanded, because - in such a case, there is no need to store the Petri net anymore. However, - in general you should assume that this field is optional, even on nodes that - are not expanded yet. - """ - return cast(nx.DiGraph, self.dag.nodes[node_id]["petri_net"]) - def edge_stable_motif( self, parent_id: int, child_id: int, reduced: bool = False ) -> BooleanSpace: """ - Return the *stable motif* associated with the specified parent-child - edge. If `reduced` is set to `False` (default), the unpercolated stable + Return the stable motif for the given parent-child edge. + + If `reduced` is set to `False` (default), the unpercolated stable motif trap space corresponding to the child node is returned; this includes the nodes that are fixed in the percolated trap space of the parent node. If `reduced` is set to `True`, the nodes fixed in the parent are removed (and thus the reduced stable motif is not a trap space of the original network, but is a maximal trap space in the network reduced by the parent node). + + Parameters + ---------- + parent_id: int + The ID of the parent node. + child_id: int + The ID of the child node. + reduced: bool + Whether to return the reduced stable motif. + + Returns + ------- + BooleanSpace + The stable motif (maximal trap space) represented by the edge. """ if reduced: @@ -472,7 +616,7 @@ def edge_stable_motif( { k: v for k, v in self.dag.edges[parent_id, child_id]["motif"].items() # type: ignore - if k not in self.node_space(parent_id) + if k not in self.node_data(parent_id)["space"] }, ) else: @@ -500,13 +644,14 @@ def expand_bfs( ) -> bool: """ Explore the succession diagram in a BFS manner. - - If `node_id` is given, initiate BFS from this node. Otherwise use - root. - - If `bfs_level_limit` is given, this is the last "level" (distance - from the initial node) - of nodes that should be expanded (any subsequent child nodes are - left unexplored). - If `size_limit` is given, the procedure stops - once `SuccessionDiagram` exceeds the given size. + + If `node_id` is given, initiate BFS from this node. Otherwise use root. + If `bfs_level_limit` is given, this is the last "level" (distance from + the initial node) of nodes that should be expanded (any subsequent child + nodes are left unexplored). + + If `size_limit` is given, the procedure stops once `SuccessionDiagram` + exceeds the given size. With default settings, the method will explore the whole succession diagram without any restrictions. @@ -615,14 +760,14 @@ def _update_node_petri_net(self, node_id: int, parent_id: int | None): such Petri net is always empty. """ - node_space = self.node_space(node_id) + node_space = self.node_data(node_id)["space"] if len(node_space) == self.network.variable_count(): # If fixed point, no need to compute. return if parent_id is not None: - parent_pn = self.node_restricted_petri_net(parent_id) + parent_pn = self.node_data(parent_id)["petri_net"] if parent_pn is None: pn = self.petri_net else: @@ -669,7 +814,9 @@ def _expand_one_node(self, node_id: int): current_space = node["space"] if DEBUG: - print(f"[{node_id}] Expanding: {len(self.node_space(node_id))} fixed vars.") + print( + f"[{node_id}] Expanding: {len(self.node_data(node_id)['space'])} fixed vars." + ) if len(current_space) == self.network.variable_count(): # This node is a fixed-point. Trappist would just @@ -686,7 +833,7 @@ def _expand_one_node(self, node_id: int): source_nodes = extract_source_variables(self.petri_net) sub_spaces: list[BooleanSpace] - pn = self.node_restricted_petri_net(node_id) + pn = self.node_data(node_id)["petri_net"] if pn is not None: # We have a pre-propagated PN for this sub-space, hence we can use # that to compute the trap spaces. @@ -747,9 +894,10 @@ def _ensure_node(self, parent_id: int | None, stable_motif: BooleanSpace) -> int child_id = None if key not in self.node_indices: child_id = self.dag.number_of_nodes() + + # Note: this must match the fields of the `NodeData` class self.dag.add_node( # type: ignore child_id, - id=child_id, # In case we ever need it within the "node data" dictionary. space=fixed_vars, depth=0, expanded=False, diff --git a/balm/types.py b/balm/types.py index 857c96b..6380c9a 100644 --- a/balm/types.py +++ b/balm/types.py @@ -13,26 +13,72 @@ class SuccessionDiagramState(TypedDict): - """A `TypedDict` class that stores the state of a succession diagram (see :class:`balm.SuccessionDiagram`). - - Parameters - ---------- - network_rules : str - The network rules as a string. - petri_net : networkx.DiGraph - The Petri net representation of the network rules. - nfvs : list[str] | None - The negative feedback vertex set used for attractor detection. - dag : networkx.DiGraph - The directed acyclic graph (DAG) representation of the succession - diagram structure. - node_indices : dict[int, int] - A dictionary mapping subspace keys to their positions in the succession - diagram (see :func:`balm.space_utils.space_unique_key`). + """ + A `TypedDict` class that stores the state of a succession diagram (see :class:`balm.SuccessionDiagram`). """ network_rules: str + """ + The network rules as a string. + """ + petri_net: nx.DiGraph + """ + The Petri net representation of the network rules. + """ + nfvs: list[str] | None + """ + The negative feedback vertex set used for attractor detection. + """ + dag: nx.DiGraph + """ + The directed acyclic graph representation of the succession diagram structure. + """ + node_indices: dict[int, int] + """ + A dictionary mapping subspace keys to their positions in the succession + diagram (see :func:`balm.space_utils.space_unique_key`). + """ + + +class NodeData(TypedDict): + """ + A `TypedDict` class that stores the data of a succession diagram node (see :class:`balm.SuccessionDiagram`). + + Returned from :func:`balm.SuccessionDiagram.node_data`. + However, this class is not directly used at runtime, and only exists for static type-checking. + Instead, at runtime, an untyped dictionary is used because that is what is returned by `networkx.DiGraph.nodes(data=True)`. + """ + + depth: int + """ + The depth of the node in the succession diagram rooted directed acyclic + graph. In cases of multiple paths from the root to the reference node, + the longest is used. + """ + + attractors: list[BooleanSpace] | None + """ + Attractor seed states for the node. If `None`, these have not been + computed, and the node may or may not have associated attractors. + """ + + petri_net: nx.DiGraph | None + """ + The Petri net representation of the node. If `None`, this have not been + computed yet. + """ + space: BooleanSpace + """ + The sub-space that the node represents (this subspace will always be a + trap space). + """ + + expanded: bool + """ + Whether the node has been expanded yet or not (`balm` builds the + succession diagram lazily). + """ diff --git a/tests/control_test.py b/tests/control_test.py index de9f4f1..37c01c3 100644 --- a/tests/control_test.py +++ b/tests/control_test.py @@ -101,7 +101,7 @@ def test_basic_succession_finding(): def test_internal_succession_control(): - sd = SuccessionDiagram.from_bnet( + sd = SuccessionDiagram.from_rules( """ S, S A, S | B @@ -143,7 +143,7 @@ def test_internal_succession_control(): def test_all_succession_control(): - sd = SuccessionDiagram.from_bnet( + sd = SuccessionDiagram.from_rules( """ S, S A, S | B @@ -185,7 +185,7 @@ def test_all_succession_control(): def test_forbidden_drivers(): - sd = SuccessionDiagram.from_bnet( + sd = SuccessionDiagram.from_rules( """ A, B & C B, A & C @@ -255,7 +255,7 @@ def test_forbidden_drivers(): def test_size_restriction(): - sd = SuccessionDiagram.from_bnet( + sd = SuccessionDiagram.from_rules( """ A, B & C B, A & C diff --git a/tests/drivers_test.py b/tests/drivers_test.py index 1e9efe6..bd5997d 100644 --- a/tests/drivers_test.py +++ b/tests/drivers_test.py @@ -1,6 +1,7 @@ -from balm.drivers import find_single_drivers, find_single_node_LDOIs from biodivine_aeon import BooleanNetwork +from balm.drivers import find_single_drivers, find_single_node_LDOIs + def test_find_single_node_LDOIs(): bn = BooleanNetwork.from_bnet( diff --git a/tests/motif_avoidant_test.py b/tests/motif_avoidant_test.py index f4c9b48..2f24c40 100644 --- a/tests/motif_avoidant_test.py +++ b/tests/motif_avoidant_test.py @@ -1,4 +1,4 @@ -from biodivine_aeon import BooleanNetwork, AsynchronousGraph +from biodivine_aeon import AsynchronousGraph, BooleanNetwork from balm.motif_avoidant import _filter_candidates # type: ignore from balm.motif_avoidant import _Pint_reachability # type: ignore diff --git a/tests/petri_net_translation_test.py b/tests/petri_net_translation_test.py index f5e3f3b..4bfa5eb 100644 --- a/tests/petri_net_translation_test.py +++ b/tests/petri_net_translation_test.py @@ -1,6 +1,6 @@ # type: ignore import pytest -from biodivine_aeon import BooleanNetwork, RegulatoryGraph +from biodivine_aeon import BooleanNetwork from networkx import DiGraph, is_isomorphic # type: ignore from balm.petri_net_translation import ( diff --git a/tests/source_SCC_test.py b/tests/source_SCC_test.py index 54f6b61..e0a0c0f 100644 --- a/tests/source_SCC_test.py +++ b/tests/source_SCC_test.py @@ -80,7 +80,7 @@ def test_find_scc_sd(): ) for node in scc_sd.node_ids(): - print(scc_sd.node_space(node)) + print(scc_sd.node_data(node)["space"]) assert scc_sd.dag.nodes[0]["space"] == {} assert scc_sd.dag.nodes[1]["space"] == {"A": 0, "B": 0} @@ -375,8 +375,8 @@ def test_isomorph(): sd_scc = SuccessionDiagram(bn) expand_source_SCCs(sd_scc) - assert [sd_bfs.node_space(id) for id in sd_bfs.node_ids()] == [ - sd_scc.node_space(id) for id in sd_scc.node_ids() + assert [sd_bfs.node_data(id)["space"] for id in sd_bfs.node_ids()] == [ + sd_scc.node_data(id)["space"] for id in sd_scc.node_ids() ] assert sd_scc.is_isomorphic(sd_bfs) diff --git a/tests/space_utils_test.py b/tests/space_utils_test.py index 739c91a..73dbf22 100644 --- a/tests/space_utils_test.py +++ b/tests/space_utils_test.py @@ -1,13 +1,13 @@ -from biodivine_aeon import BooleanNetwork, BooleanExpression, AsynchronousGraph +from biodivine_aeon import AsynchronousGraph, BooleanExpression, BooleanNetwork from balm.space_utils import ( expression_to_space_list, is_subspace, - percolate_network, - percolation_conflicts, percolate_expression, + percolate_network, percolate_space, percolate_space_strict, + percolation_conflicts, space_unique_key, ) diff --git a/tests/succession_diagram_test.py b/tests/succession_diagram_test.py index 577526b..9418cb0 100644 --- a/tests/succession_diagram_test.py +++ b/tests/succession_diagram_test.py @@ -30,7 +30,7 @@ def test_succession_diagram_structure(self): assert ( max( [ - succession_diagram.node_depth(i) + succession_diagram.node_data(i)["depth"] for i in succession_diagram.node_ids() ] ) @@ -80,7 +80,7 @@ def test_succession_diagram_structure(self): assert ( max( [ - succession_diagram.node_depth(i) + succession_diagram.node_data(i)["depth"] for i in succession_diagram.node_ids() ] ) @@ -109,8 +109,8 @@ def test_succession_diagram_structure(self): def test_state(): - sd1 = SuccessionDiagram.from_bnet("A, B\nB, A") - sd2 = SuccessionDiagram.from_bnet("A, B\nB, A\nC, C & B") + sd1 = SuccessionDiagram.from_rules("A, B\nB, A") + sd2 = SuccessionDiagram.from_rules("A, B\nB, A\nC, C & B") sd1.build() sd2.__setstate__(sd1.__getstate__()) slots1 = [x + str(sd1.__getattribute__(x)) for x in sd1.__slots__] @@ -193,7 +193,7 @@ def test_expansion_comparisons(network_file: str): # This should always create a succession_diagram with exactly one minimal trap space, # as the rest for min_trap in sd_bfs.minimal_trap_spaces(): - space = sd_bfs.node_space(min_trap) + space = sd_bfs.node_data(min_trap)["space"] sd_target = SuccessionDiagram(bn) assert sd_target.expand_to_target(space, size_limit=NODE_LIMIT) @@ -308,3 +308,16 @@ def test_attractor_expansion(network_file: str): # All symbolic attractors must be covered by some seed at this point. assert len(symbolic_attractors) == 0 + + +def test_attractor_extraction(): + sd = balm.SuccessionDiagram.from_rules( + """ + A, B + B, A & C + C, !A | B + """ + ) + sd.build() + eas = sd.expanded_attractor_seeds() + assert eas == {1: [{"A": 0, "B": 0, "C": 1}], 2: [{"A": 1, "B": 1, "C": 1}]} diff --git a/tests/trappist_core_test.py b/tests/trappist_core_test.py index 33eda6c..cb77e7e 100644 --- a/tests/trappist_core_test.py +++ b/tests/trappist_core_test.py @@ -1,12 +1,8 @@ -from biodivine_aeon import ( - BooleanNetwork, - FixedPoints, - AsynchronousGraph, -) +from biodivine_aeon import AsynchronousGraph, BooleanNetwork, FixedPoints +from balm.interaction_graph_utils import cleanup_network from balm.petri_net_translation import network_to_petrinet from balm.trappist_core import compute_fixed_point_reduced_STG, trappist -from balm.interaction_graph_utils import cleanup_network from balm.types import BooleanSpace