Boolean Network

This module defines the BooleanNetwork class, which provides a high-level framework for modeling, simulating, and analyzing Boolean networks.

A BooleanNetwork represents a discrete dynamical system \(F = (f_1, \ldots, f_N)\) composed of multiple BooleanFunction objects as update rules. The class includes methods for constructing state transition graphs, identifying attractors, computing robustness and sensitivity measures, and exporting truth tables.

Several computational routines—particularly those involving state space exploration, attractor detection, and robustness estimation—offer optional Numba-based just-in-time (JIT) acceleration. Installing Numba is recommended for optimal performance but not required; all features remain functional without it.

This module serves as the central interface for dynamic Boolean network analysis within the BoolForge package.

Example

>>> from boolforge import BooleanNetwork
>>> bn = BooleanNetwork(F=[[0, 1], [0, 0, 0, 1], [0, 1]], I=[[1], [0, 2], [1]])
>>> bn.get_attractors_synchronous_exact()
boolforge.boolean_network.get_entropy_of_basin_size_distribution(basin_sizes: Sequence[int] | ndarray) → float[source]

Compute the Shannon entropy of a basin size distribution.

The basin sizes are first normalized to form a probability distribution. The Shannon entropy is then computed as

H = -sum(p_i * log(p_i)),

where p_i is the proportion of states in basin i.

Parameters

basin_sizesSequence[int] or np.ndarray

Sizes of the basins of attraction, where each entry corresponds to the number of initial conditions that converge to a given attractor.

Returns

float

Shannon entropy of the basin size distribution.

class boolforge.boolean_network.BooleanNetwork(F: Sequence[BooleanFunction | list[int] | ndarray], I: Sequence[Sequence[int]] | WiringDiagram, variables: Sequence[str] | None = None, SIMPLIFY_FUNCTIONS: bool = False)[source]

Bases: WiringDiagram

Representation of a Boolean network.

A Boolean network consists of a wiring diagram specifying regulatory interactions between nodes and a collection of Boolean update functions defining the dynamics at each node.

In a BooleanNetwork, constant nodes are removed during initialization, so all nodes represent dynamic variables.

Parameters

Fsequence

Sequence of Boolean update functions or truth tables. Each entry may be a BooleanFunction instance, a truth table, or a Boolean expression. The length of F must match the number of nodes in the wiring diagram.

Isequence of sequences of int or WiringDiagram

Wiring diagram specifying the regulators of each node, or an existing WiringDiagram instance.

variablessequence of str, optional

Names of the variables corresponding to each node. Ignored if I is provided as a WiringDiagram.

SIMPLIFY_FUNCTIONSbool, optional

If True, simplify Boolean update functions after initialization. Default is False.

Attributes

Flist[BooleanFunction]

Boolean update functions for each node.

Ilist[np.ndarray[int]]

Wiring diagram specifying the regulators of each node.

variablesnp.ndarray[str]

Names of the variables corresponding to each node.

Nint

Number of dynamic (non-constant) nodes in the network.

indegreesnp.ndarray[int]

Indegree of each node.

outdegreesnp.ndarray[int]

Outdegree of each node.

constantsdict[str, dict[str, int | list[str]]]

Mapping of node indices to constant values.

weightslist[np.ndarray[float]] or None

Interaction weights associated with the wiring diagram.

STGdict or None

State transition graph, initialized to None and computed on demand.

remove_constants() → None[source]

Remove structurally constant nodes from the Boolean network.

A node is considered constant if it has no regulators (indegree zero). Such nodes are eliminated from the dynamic network by propagating their fixed Boolean values to downstream nodes. Eliminated constants and their effects are recorded in the constants attribute.

Notes

  • The Boolean value of a constant node is taken from its Boolean function.

  • After removal, self.N refers to the number of remaining dynamic nodes.

  • Nodes that lose all regulators as a result of constant removal are assigned a non-essential self-loop to preserve network structure.

classmethod from_cana(cana_BooleanNetwork: BooleanNetwork) → BooleanNetwork[source]

Construct a BooleanNetwork from a cana.BooleanNetwork instance.

This compatibility method converts a Boolean network defined using the cana package into a BoolForge BooleanNetwork.

Parameters

cana_BooleanNetworkcana.boolean_network.BooleanNetwork

A Boolean network instance from the cana package.

Returns

BooleanNetwork

The corresponding BoolForge BooleanNetwork.

Raises

TypeError

If the input object does not appear to be a valid CANA BooleanNetwork.

KeyError

If required fields are missing from the CANA logic specification.

classmethod from_string(network_string: str, separator: str | Sequence[str] = ',', max_degree: int = 24, original_not: str | Sequence[str] = 'NOT', original_and: str | Sequence[str] = 'AND', original_or: str | Sequence[str] = 'OR', ALLOW_TRUNCATION: bool = False) → BooleanNetwork[source]

Construct a BooleanNetwork from a textual Boolean rule specification.

This compatibility method parses a string representation of Boolean update rules and constructs a corresponding BooleanNetwork. The input format is intended for legacy or trusted sources and supports logical expressions using AND/OR/NOT operators.

Warning

This method uses eval internally and MUST NOT be used on untrusted input. It is provided solely for backward compatibility.

Parameters

network_stringstr

String encoding Boolean update rules, one per line.

separatorstr or sequence of str, optional

Separator(s) between variable names and Boolean expressions.

max_degreeint, optional

Maximum allowed indegree for explicit truth-table construction.

original_not, original_and, original_orstr or sequence of str, optional

Operator strings to be replaced by logical NOT, AND, OR.

ALLOW_TRUNCATIONbool, optional

If False (default), nodes with indegree greater than max_degree raise a ValueError. If True, such nodes are replaced by identity self-loops, allowing fast construction of large networks while ignoring high-degree functions.

Returns

BooleanNetwork

The constructed Boolean network.

Raises

ValueError

If parsing fails or if ALLOW_TRUNCATION is False and a node exceeds max_degree.

classmethod from_DiGraph(nx_DiGraph: DiGraph) → WiringDiagram[source]

Construct a WiringDiagram from a NetworkX directed graph.

Each node in the directed graph represents a Boolean variable, and each directed edge u -> v indicates that variable u regulates variable v.

Parameters

nx_DiGraphnx.DiGraph

Directed graph whose nodes represent variables and whose edges represent regulatory interactions.

Node attributes (optional)
namestr

Name of the variable. If not provided, the node label is used when possible.

Edge attributes (optional)
weightfloat

Weight of the regulatory interaction from u to v. If present on all edges, weights are stored in the weights attribute of the resulting WiringDiagram.

Returns

WiringDiagram

Wiring diagram representing the topology of the directed graph.

Notes

  • Nodes are ordered according to the iteration order of nx_DiGraph.nodes.

  • Regulator lists are constructed from incoming edges (graph predecessors).

  • Edge weights are only stored if all edges define a 'weight' attribute.

Examples

>>> import networkx as nx
>>> from boolforge import WiringDiagram
>>> G = nx.DiGraph()
>>> G.add_edges_from([(0, 1), (1, 2)])
>>> W = WiringDiagram.from_DiGraph(G)
>>> W.I
[array([], dtype=int64), array([0]), array([1])]
>>> W.get_source_nodes()
{0: True, 1: False, 2: False}
to_cana() → BooleanNetwork[source]

Export the Boolean network as a cana.BooleanNetwork instance.

This compatibility method converts the current BooleanNetwork into an equivalent representation from the cana package. The exported network reflects the current state of the model, including any removed constants, simplifications, or identity self-loops.

Returns

cana.boolean_network.BooleanNetwork

A cana BooleanNetwork instance representing this network.

Raises

ImportError

If the cana package is not installed.

to_bnet(separator: str = ',\t', AS_POLYNOMIAL: bool = True) → str[source]

Export the Boolean network in BNET format.

This compatibility method returns a string representation of the Boolean network in the BNET format used by tools such as BoolNet and PyBoolNet, with one line per variable of the form variable <separator> function.

Parameters

separatorstr, optional

String used to separate the target variable from its update function. Default is “,t”.

AS_POLYNOMIALbool, optional

If True (default), return Boolean functions in polynomial form. If False, return functions as logical expressions.

Returns

str

A string containing the BNET representation of the network.

to_truth_table(filename: str | None = None) → DataFrame[source]

Construct the full synchronous truth table of the Boolean network.

Each row corresponds to a network state at time t and its deterministic successor at time t+1 under synchronous updating.

Parameters

filenamestr, optional

If provided, the truth table is written to a file. The file extension determines the format and must be one of 'csv', 'xls', or 'xlsx'. If None (default), no file is created.

Returns

pandas.DataFrame

The full truth table with shape (2**N, 2*N).

Notes

  • States are enumerated in lexicographic order, consistent with utils.get_left_side_of_truth_table.

  • This method computes and stores the synchronous state transition graph (self.STG) if it has not been computed previously.

  • Exporting to Excel requires the openpyxl package.

__str__()[source]

Return str(self).

__repr__()[source]

Return repr(self).

__call__(state)[source]

Apply one synchronous update step to the Boolean network.

The next state is obtained by evaluating each node’s Boolean update function on the current values of its regulators.

Parameters

statesequence of int

Current network state as a binary vector of length N, ordered according to self.variables.

Returns

np.ndarray

The updated network state after one synchronous update.

Notes

This method is equivalent to calling update_network_synchronously.

get_types_of_regulation() → list[ndarray][source]

Compute and return regulation types (weights) for all nodes in the network.

For each Boolean function, the type of each input regulation is determined via BooleanFunction.get_type_of_inputs and mapped to numerical weights using dict_weights. The resulting weights are stored in the self.weights attribute and also returned.

Returns

list of np.ndarray

Regulation weights for each node, aligned with the wiring diagram.

Notes

  • This method recomputes self.weights from scratch.

  • Calling this method overwrites any existing values in self.weights.

simplify_functions() → None[source]

Remove all non-essential regulators from the Boolean network.

For each node, non-essential regulators (identified via np.nan entries in self.weights) are removed from the wiring diagram and the associated Boolean function is restricted to its essential inputs. Nodes that would otherwise lose all regulators are assigned an identity self-loop to preserve network structure.

Notes

  • This method modifies the network in place.

  • Regulation types (self.weights) are recomputed if necessary.

  • Identity self-loops introduced here are structural artifacts and do not represent genuine regulatory interactions.

get_identity_nodes(AS_DICT: bool = False) → dict[int, bool] | ndarray[source]

Identify identity (memory) nodes in the Boolean network.

An identity node is a node with a single self-regulatory edge whose Boolean update function is the identity function f(x) = x. Such nodes retain their state over time unless externally modified.

Parameters

AS_DICTbool, optional

If True, return a dictionary mapping node indices to booleans. If False (default), return an array of indices of identity nodes.

Returns

dict[int, bool] or np.ndarray

If AS_DICT is True, a dictionary indicating which nodes are identity nodes. If AS_DICT is False, an array of indices of identity nodes.

get_network_with_fixed_identity_nodes(values_identity_nodes: Sequence[int]) → BooleanNetwork[source]

Construct a Boolean network with identity nodes fixed to given values.

Identity nodes are nodes with a single self-regulatory edge and identity update rule f(x) = x. This method fixes the values of such nodes and returns a new BooleanNetwork with the corresponding constants removed.

Parameters

values_identity_nodessequence of int

Values to fix for each identity node, in the order returned by get_identity_nodes(AS_DICT=False). Each value must be either 0 or 1.

Returns

BooleanNetwork

A new BooleanNetwork with the specified identity nodes fixed. Any constants previously removed from the original network are preserved.

get_network_with_node_controls(indices_controlled_nodes: Sequence[int], values_controlled_nodes: Sequence[int], KEEP_CONTROLLED_NODES: bool = False) → BooleanNetwork[source]

Construct a Boolean network with specified nodes fixed to given values.

This method applies node-level interventions by fixing selected nodes to constant Boolean values. Controlled nodes may either be removed from the dynamic network as constants or retained as identity-clamped nodes.

Parameters

indices_controlled_nodessequence of int

Indices of nodes to be fixed.

values_controlled_nodessequence of int

Values to fix for each specified node, in the same order as indices_controlled_nodes. Each value must be either 0 or 1.

KEEP_CONTROLLED_NODESbool, optional

If True, controlled nodes are retained in the network as identity nodes with self-loops. If False (default), controlled nodes are eliminated as constants.

Returns

BooleanNetwork

A new BooleanNetwork with the specified node controls applied.

get_network_with_edge_controls(control_targets: Sequence[int], control_sources: Sequence[int], values_edge_controls: Sequence[int] | None = None) → BooleanNetwork[source]

Construct a Boolean network with specified regulatory edges controlled.

This method fixes the influence of selected source nodes on selected target nodes by restricting the target’s Boolean update function to entries where the source assumes a specified value, and then removing the corresponding regulatory edge.

Parameters

control_targetssequence of int

Indices of target nodes.

control_sourcessequence of int

Indices of source nodes whose influence on the corresponding targets is to be controlled.

values_edge_controlssequence of int, optional

Fixed values (0 or 1) imposed on each controlled edge. If None, all controlled edges are fixed to 0.

Returns

BooleanNetwork

A new BooleanNetwork with the specified edge controls applied.

Raises

ValueError

If input lengths do not match, indices are invalid, or edge values are not in {0, 1}.

update_single_node(index: int, states_regulators: Sequence[int]) → int[source]

Update the state of a single node.

The new state is obtained by applying the Boolean update function to the states of its regulators.

Parameters

indexint

Index of the node to update.

states_regulatorssequence of int

Binary states of the node’s regulators.

Returns

int

Updated state of the node (0 or 1).

update_network_synchronously(state: Sequence[int]) → ndarray[source]

Perform a synchronous update of the Boolean network.

Parameters

statesequence of int

Binary state vector of length N.

Returns

np.ndarray

Updated state vector.

update_network_SDDS(state: Sequence[int], P: ndarray, *, rng=None) → ndarray[source]

Perform a stochastic discrete dynamical system (SDDS) update of the network.

This update scheme follows the SDDS formalism: for each node, the deterministic Boolean update is first computed. If the update would increase the node’s state, the change occurs with the node-specific activation probability. If the update would decrease the node’s state, the change occurs with the node-specific degradation probability. Otherwise, the node’s state remains unchanged.

Parameters

statesequence of int

Current network state (binary vector of length N).

Pnp.ndarray

Array of shape (N, 2), where P[i, 0] is the activation probability and P[i, 1] is the degradation probability for node i.

rngoptional

Random number generator or seed, passed to utils._coerce_rng.

Returns

np.ndarray

Updated network state after one stochastic SDDS update.

Notes

This implementation follows the SDDS framework introduced in:

Murrugarra, D., Veliz-Cuba, A., Aguilar, B., Arat, S., & Laubenbacher, R. (2012). Modeling stochasticity and variability in gene regulatory networks. EURASIP Journal on Bioinformatics and Systems Biology, 2012(1), 5.

The method assumes that state is a valid binary vector and that P has the correct shape; no additional validation is performed for performance reasons.

get_steady_states_asynchronous_exact(stochastic_weights: Sequence[float] | None = None, max_iterations: int = 1000, tol: float = 1e-09) → dict[source]

Compute the steady states and basin probabilities under general asynchronous update.

This method exhaustively constructs the asynchronous state transition graph (STG) of the Boolean network under a general asynchronous update scheme, where nodes are selected for update according to given propensities. The resulting Markov chain is solved exactly using an iterative Gauss–Seidel scheme to obtain absorption probabilities into steady states.

Parameters

stochastic_weightssequence of float or None, optional

Relative update propensities for each node. If None (default), all nodes are updated with equal probability. The weights are normalized internally.

max_iterationsint, optional

Maximum number of Gauss–Seidel iterations used to compute absorption probabilities before declaring non-convergence.

tolfloat, optional

Convergence tolerance for the infinity norm of probability updates.

s

Returns

dict

Dictionary with the following entries:

  • "SteadyStates" : list of int Decimal representations of steady states.

  • "NumberOfSteadyStates" : int Total number of steady states.

  • "BasinSizes" : np.ndarray Fraction of the state space converging to each steady state.

  • "STGAsynchronous" : dict Asynchronous state transition graph represented as a Markov kernel.

  • "FinalTransitionProbabilities" : np.ndarray Absorption probabilities from each state into each steady state.

Raises

RuntimeError

If the iterative solver does not converge within max_iterations.

get_steady_states_asynchronous(nsim: int = 500, initial_states: Sequence[int] | None = None, search_depth: int = 50, DEBUG: bool = False, *, rng=None) → dict[source]

Approximate steady states of a Boolean network under asynchronous updates.

This method performs a Monte Carlo–style exploration of the asynchronous state space by simulating asynchronous updates from a collection of initial states. Each simulation proceeds until a steady state is reached or until a maximum search depth is exceeded.

Unlike get_steady_states_asynchronous_exact, this method does not exhaustively explore the full state space and does not guarantee that all steady states will be found. It is intended for large networks where exact enumeration is infeasible.

Parameters

nsimint, optional

Number of asynchronous simulations to perform (default is 500).

initial_statessequence of int or None, optional

Initial states to use for the simulations, given as decimal representations of network states. If None (default), nsim random initial states are generated.

search_depthint, optional

Maximum number of asynchronous update steps per simulation before giving up on convergence (default is 50).

DEBUGbool, optional

If True, print detailed debugging information during simulation.

rngoptional

Random number generator or seed, passed to utils._coerce_rng.

Returns

dict

Dictionary with the following entries:

  • "SteadyStates" : list of int Decimal representations of steady states encountered.

  • "NumberOfSteadyStates" : int Number of unique steady states found.

  • "BasinSizes" : list of int Counts of how many simulations converged to each steady state.

  • "STGAsynchronous" : dict Partial cache of asynchronous transitions encountered during simulation. Keys are (state, node_index) and values are successor states (all in decimal form).

  • "InitialSamplePoints" : list of int Decimal initial states used in the simulations (either provided explicitly or generated randomly).

Notes

  • This method detects only steady states (fixed points). If the asynchronous dynamics contain limit cycles, simulations may fail to converge within search_depth.

  • The returned asynchronous transition graph is generally incomplete and should be interpreted as a cache of explored transitions rather than the full STG.

get_steady_states_asynchronous_given_one_initial_condition(initial_condition: int | Sequence[int] = 0, nsim: int = 500, stochastic_weights: Sequence[float] | None = None, search_depth: int = 50, DEBUG: bool = False, *, rng=None) → dict[source]

Approximate steady states reachable from a single initial condition under asynchronous updates.

This method performs multiple asynchronous simulations starting from the same initial condition. In each simulation, nodes are updated one at a time according to either a uniform random order or node-specific stochastic update propensities. The simulation proceeds until a steady state is reached or a maximum number of update steps is exceeded.

The method is sampling-based and does not guarantee that all reachable steady states are found. It is intended for exploratory analysis and for networks where exhaustive asynchronous analysis is infeasible.

Parameters

initial_conditionint or sequence of int, optional

Initial network state. If an integer is provided, it is interpreted as the decimal encoding of a Boolean state. If a sequence is provided, it must be a binary vector of length N. Default is 0.

nsimint, optional

Number of asynchronous simulation runs (default is 500).

stochastic_weightssequence of float or None, optional

Relative update propensities for each node. If provided, must have length N and be strictly positive. The weights are normalized internally. If None (default), nodes are updated uniformly at random.

search_depthint, optional

Maximum number of asynchronous update steps per simulation.

DEBUGbool, optional

If True, print detailed debugging information during simulation.

rngoptional

Random number generator or seed, passed to utils._coerce_rng.

Returns

dict

Dictionary with the following entries:

  • "SteadyStates" : list of int Decimal representations of steady states reached.

  • "NumberOfSteadyStates" : int Number of unique steady states found.

  • "BasinSizes" : list of int Number of simulations converging to each steady state.

  • "TransientTimes" : list of list of int For each steady state, a list of transient lengths (number of update steps before convergence).

  • "STGAsynchronous" : dict Partial cache of asynchronous transitions encountered during simulation. Keys are (state, node_index) and values are successor states (all in decimal form).

  • "UpdateQueues" : list of list of int For each simulation, the sequence of visited states (in decimal form).

Notes

  • Only steady states (fixed points) are detected. If the asynchronous dynamics contain limit cycles, simulations may fail to converge within search_depth.

  • The returned asynchronous transition graph is incomplete and represents only transitions encountered during sampling.

get_attractors_synchronous(nsim: int = 500, initial_sample_points: Sequence[int | Sequence[int]] | None = None, n_steps_timeout: int = 1000, INITIAL_SAMPLE_POINTS_AS_BINARY_VECTORS: bool = False, USE_NUMBA: bool = True, *, rng=None) → dict[source]

Approximate synchronous attractors of a Boolean network via sampling.

This method estimates the synchronous attractors (fixed points and cycles) of a Boolean network by simulating synchronous updates from a collection of initial states. For each simulation, the network is updated until an attractor is reached or a maximum number of update steps is exceeded.

The method is sampling-based and does not guarantee that all attractors are found. Basin sizes are lower-bound estimates based on the sampled initial conditions.

If Numba is available and USE_NUMBA=True, synchronous updates are accelerated using a compiled kernel.

Parameters

nsimint, optional

Number of initial conditions to simulate (default is 500). Ignored if initial_sample_points is provided.

initial_sample_pointssequence of int or sequence of sequence of int, optional

Initial states to use. If provided, its length determines the number of simulations. Interpretation depends on INITIAL_SAMPLE_POINTS_AS_BINARY_VECTORS.

n_steps_timeoutint, optional

Maximum number of synchronous update steps per simulation before declaring a timeout (default is 1000).

INITIAL_SAMPLE_POINTS_AS_BINARY_VECTORSbool, optional

If True, initial_sample_points are interpreted as binary vectors; otherwise they are interpreted as decimal-encoded states.

USE_NUMBAbool, optional

If True (default) and Numba is available, use a Numba-accelerated synchronous update kernel.

rngoptional

Random number generator or seed, passed to utils._coerce_rng.

Returns

dict

Dictionary with the following entries:

  • "Attractors" : list of list of int Attractors found, each represented as a list of decimal states (cycles are given in cyclic order).

  • "NumberOfAttractors" : int Number of distinct attractors found.

  • "BasinSizes" : list of int Number of sampled initial conditions converging to each attractor.

  • "AttractorDict" : dict Mapping from visited states (decimal) to attractor index.

  • "InitialSamplePoints" : list of int Decimal initial states used for sampling.

  • "STG" : dict Sampled synchronous state transition graph (state → successor state).

  • "NumberOfTimeouts" : int Number of simulations that did not converge within n_steps_timeout.

Notes

  • This method is intended for networks with long transient dynamics, where exhaustive synchronous analysis is infeasible.

  • Basin sizes are sampling-based estimates and should not be interpreted as exact proportions of the state space.

compute_synchronous_state_transition_graph(USE_NUMBA: bool = True) → None[source]

Compute the exact synchronous state transition graph (STG).

The STG is stored in self.STG as a one-dimensional NumPy array of length 2**N, where self.STG[x] is the decimal representation of the successor state reached from state x under synchronous updating.

This computation is exact and requires memory proportional to 2**N. It is therefore intended for small-to-moderate networks only.

Parameters

USE_NUMBAbool, optional

If True (default) and Numba is available, use a compiled kernel to accelerate computation.

get_attractors_synchronous_exact(USE_NUMBA: bool = True) → dict[source]

Compute all attractors and their exact basin sizes under synchronous updating.

This method computes the exact synchronous state transition graph (STG) and analyzes it as a functional graph on 2**N states. All attractors (cycles), their basin sizes, and the attractor reached from each state are determined exactly.

This computation requires memory and time proportional to 2**N and is intended for small-to-moderate networks only.

Parameters

USE_NUMBAbool, optional

If True (default) and Numba is available, use a compiled kernel for attractor detection.

Returns

dict

Dictionary with keys:

  • Attractorslist[list[int]]

    Each attractor represented as a list of decimal states forming a cycle.

  • NumberOfAttractorsint

    Total number of attractors.

  • BasinSizesnp.ndarray[float]

    Fraction of all states belonging to each attractor basin.

  • AttractorIDnp.ndarray[int]

    For each of the 2**N states, the index of the attractor it reaches.

  • STGnp.ndarray[int]

    The synchronous state transition graph.

get_transient_lengths_exact(USE_NUMBA: bool = True) → ndarray[source]
Compute exact transient length using:

  • Full STG from get_attractors_synchronous_exact()

  • Attractors (cycle states) from get_attractors_synchronous_exact()

This avoids indegree-pruning because cycle states are given explicitly.

get_attractors_and_robustness_measures_synchronous_exact(USE_NUMBA: bool = True, GET_STRATIFIED_COHERENCES: bool = False) → dict[source]

Compute attractors and exact robustness measures of a synchronously updated Boolean network.

This method constructs the exact synchronous state transition graph (STG) on 2**N states and analyzes it as a functional graph. All attractors (cycles), basin sizes, and the attractor reached from each state are determined exactly. Based on this decomposition, exact coherence and fragility measures are computed for the full network, for each basin of attraction, and for each attractor.

Optionally, coherence can be stratified by the transient length (distance from the attractor) of each state, allowing robustness to be analyzed as a function of how far states lie from their eventual attractor.

This computation requires memory and time proportional to 2**N and is intended for small-to-moderate networks. When Numba is enabled, exact and stratified robustness measures remain feasible up to moderate values of N (e.g., N ≈ 20 on typical hardware).

Parameters

USE_NUMBAbool, optional

If True (default) and Numba is available, compiled kernels are used for robustness and transient-length computations, resulting in substantial speedups.

GET_STRATIFIED_COHERENCESbool, optional

If True, coherence is additionally computed as a function of the transient length (distance to the attractor) of each state. When Numba is enabled, this option incurs only modest additional computational cost. Default is False.

Returns

dict

Dictionary with the following keys:

  • Attractorslist[list[int]]

    Each attractor represented as a list of decimal states forming a cycle.

  • NumberOfAttractorsint

    Total number of attractors.

  • BasinSizesnp.ndarray of float

    Fraction of all states belonging to each attractor basin.

  • AttractorIDnp.ndarray of int

    For each of the 2**N states, the index of the attractor it eventually reaches.

  • Coherencefloat

    Exact global network coherence.

  • Fragilityfloat

    Exact global network fragility.

  • BasinCoherencenp.ndarray of float

    Exact coherence of each basin of attraction.

  • BasinFragilitynp.ndarray of float

    Exact fragility of each basin of attraction.

  • AttractorCoherencenp.ndarray of float

    Exact coherence of each attractor.

  • AttractorFragilitynp.ndarray of float

    Exact fragility of each attractor.

If GET_STRATIFIED_COHERENCES is True, the dictionary additionally contains:

  • StratifiedCoherencesnp.ndarray of float

    Coherence values stratified by attractor and transient length.

  • DistanceFromAttractorCountnp.ndarray of int

    Number of state–hypercube-edge incidences contributing to each stratified coherence entry.

  • DistanceFromAttractornp.ndarray of int

    Transient length (distance to attractor) for each state.

get_attractors_and_robustness_measures_synchronous(number_different_IC: int = 500, RETURN_ATTRACTOR_COHERENCE: bool = True, *, rng=None) → dict[source]

Approximate attractors and robustness measures under synchronous updating.

This method samples the attractor landscape by simulating the network from multiple random initial conditions (ICs) and their single-bit perturbations. It returns Monte-Carlo approximations of global coherence, fragility, and a final Hamming-distance-based measure, along with per-basin approximations. Optionally, it additionally estimates attractor-level coherence and fragility by perturbing attractor states found during sampling.

Notes

  • The attractor set returned is a lower bound on the true number of attractors, because only the sampled portion of state space is explored.

  • For N >= 64, decimal encoding of states may exceed np.int64 and this method uses bitstrings (type str) as state identifiers.

Parameters

number_different_ICint, optional

Number of random initial conditions to sample (default is 500). For each IC, the method also simulates one randomly chosen single-bit perturbation.

RETURN_ATTRACTOR_COHERENCEbool, optional

If True (default), also compute attractor-level coherence and fragility by perturbing attractor states found during sampling.

rngNone or numpy.random.Generator, optional

Random number generator or seed-like object. Passed to utils._coerce_rng.

Returns

dict

Dictionary with keys:

  • Attractorslist[list[int]] or list[list[str]]

    List of discovered attractors, each represented as a list of states forming a cycle. States are decimals (int) for N < 64 and bitstrings (str) for N >= 64.

  • LowerBoundOfNumberOfAttractorsint

    Number of distinct attractors discovered (a lower bound on the true number of attractors).

  • BasinSizesApproximationnp.ndarray[float]

    Approximate basin size (fraction of sampled trajectories that end in each attractor). Sums to ~1 over discovered attractors.

  • CoherenceApproximationfloat

    Approximate global coherence: probability that a random IC and its single-bit perturbation reach the same attractor.

  • FragilityApproximationfloat

    Approximate global fragility: expected normalized difference between reached attractors when the IC and perturbation reach different attractors. Normalized by N.

  • FinalHammingDistanceApproximationfloat

    Approximate final Hamming distance between the two periodic trajectories when comparing the IC and its perturbation. This is a distance in [0, 1], where 0 means identical and 1 means completely different.

  • BasinCoherenceApproximationnp.ndarray[float]

    Approximate coherence per basin (same definition as coherence but conditioned on having reached that basin).

  • BasinFragilityApproximationnp.ndarray[float]

    Approximate fragility per basin (same definition as fragility but conditioned on having reached that basin).

  • AttractorCoherencenp.ndarray[float], optional

    If RETURN_ATTRACTOR_COHERENCE is True: estimated attractor-level coherence (probability that a single-bit perturbation of an attractor state returns to the same attractor).

  • AttractorFragilitynp.ndarray[float], optional

    If RETURN_ATTRACTOR_COHERENCE is True: estimated attractor-level fragility based on differences between the original attractor and the attractor reached after perturbation.

References

Park, K. H., Costa, F. X., Rocha, L. M., Albert, R., & Rozum, J. C. (2023). Models of cell processes are far from the edge of chaos. PRX Life, 1(2), 023009.

Bavisetty, V. S. N., Wheeler, M., & Kadelka, C. (2025). Attractors are less stable than their basins: Canalization creates a coherence gap in gene regulatory networks. bioRxiv 2025-11.

get_derrida_value(nsim: int = 1000, EXACT: bool = False, USE_NUMBA: bool = True, *, rng=None) → float[source]

Compute the Derrida value of a Boolean network.

The Derrida value measures the average Hamming distance between the one-step synchronous updates of two states that differ by a single-bit perturbation. It quantifies the short-term sensitivity of the network dynamics to small perturbations.

If EXACT is True, the Derrida value is computed exactly as the mean (unnormalized) average sensitivity of the Boolean update functions. Otherwise, it is approximated via Monte Carlo simulation.

Parameters

nsimint, optional

Number of Monte Carlo simulations to perform (default is 1000). Ignored if EXACT is True.

EXACTbool, optional

If True, compute the exact Derrida value. If False (default), approximate the Derrida value using Monte Carlo simulation.

USE_NUMBAbool, optional

If True (default) and Numba is available, use a compiled kernel for Monte Carlo simulation.

rngNone or np.random.Generator, optional

Random number generator, passed through utils._coerce_rng.

Returns

float

The Derrida value, defined as the average Hamming distance after one synchronous update following a single-bit perturbation.

References

Derrida, B., & Pomeau, Y. (1986). Random networks of automata: a simple annealed approximation. Europhysics Letters, 1(2), 45.