API Reference

batss: A Python package with Rust backend for simulation of chemical reaction networks through a fast batchign algorithm.

This package provides tools for simulating population protocols and chemical reaction networks. The core simulation engine is implemented in Rust for performance, with a Python interface for ease of use and visualization.

Modules:
  • batss.simulation: Module for simulating chemical reaction networks

  • batss.crn: Module for expressing population protocols using chemical reaction network notation

  • batss.snapshot: Module for visualizing species counts during or after simulation

Example

from batss import species, Simulation

# Define a simple protocol
a, b, u = species('A B U')
approx_majority = [
    a + b >> 2 * u,
    a + u >> 2 * a,
    b + u >> 2 * b,
]

# Initialize and run simulation
n = 10 ** 5
a_init = int(0.51 * n)
b_init = n - a
init_config = {a: a_init, b: b_init}
sim = Simulation(init_config=init_config, rule=approx_majority)
sim.run()
sim.history.plot()

Module Documentation

A module for simulating population protocols.

The main class Simulation is created with a description of the protocol and the initial condition, and is responsible for running the simulation.

The general syntax is

a, b, u = 'A', 'B', 'U'
approx_majority = {
    (a,b): (u,u),
    (a,u): (a,a),
    (b,u): (b,b),
}
n = 10 ** 5
init_config = {a: 0.51 * n, b: 0.49 * n}
sim = Simulation(init_config=init_config, rule=approx_majority)
sim.run()
sim.history.plot()

More examples given in https://github.com/UC-Davis-molecular-computing/ppsim/tree/main/examples

time_trials() is a convenience function used for gathering data about the convergence time of a protocol.

batss.simulation.Rule

Type representing transition rules for protocol. Is one of three types:

  • a function that takes two states and returns either a tuple of two states or a dictionary mapping pairs of states to probabilities.

  • a dictionary mapping pairs of states to either a tuple of two states or a dictionary mapping pairs of states to probabilities.

  • a list of ppsim.crn.Reaction objects describing a CRN, which will be passed as a CRN if the simulator method is CRN,

or will be converted into a population protocol if it’s using original ppsim-style multibatching.

alias of Callable[[Hashable, Hashable], tuple[Hashable, Hashable] | dict[tuple[Hashable, Hashable], float]] | dict[tuple[Hashable, Hashable], tuple[Hashable, Hashable] | dict[tuple[Hashable, Hashable], float]] | Iterable[Reaction]

batss.simulation.state_enumeration(init_dist, rule)[source]

Finds all reachable states by breadth-first search.

Parameters:
Return type:

set[Hashable]

Returns:

a set of all reachable states

class batss.simulation.Simulation(init_config, rule, *, simulator_method='crn', transition_order='symmetric', seed=None, volume=None, continuous_time=False, time_units=None, crn=None, **kwargs)[source]

Bases: object

Class to simulate a population protocol.

Parameters:

Initialize a Simulation.

Parameters:
  • init_config (dict[TypeVar(StateTypeVar, bound= Hashable), int]) – The starting configuration, as a dictionary mapping states to counts.

  • rule (Union[Callable[[Hashable, Hashable], tuple[Hashable, Hashable] | dict[tuple[Hashable, Hashable], float]], dict[tuple[Hashable, Hashable], tuple[Hashable, Hashable] | dict[tuple[Hashable, Hashable], float]], Iterable[Reaction]]) – A representation of the transition rule. The first two options are a dictionary, whose keys are tuples of 2 states and values are their outputs, or a function which takes pairs of states as input. For a deterministic transition function, the output is a tuple of 2 states. For a probabilistic transition function, the output is a dictionary mapping tuples of states to probabilities. Inputs that are not present in the dictionary, or return None from the function, are interpreted as null transitions that return the same pair of states as the output. The third option is a list of ppsim.crn.Reaction objects describing a CRN, which will be parsed into an equivalent population protocol.

  • simulator_method (str) –

    Which Simulator method to use, either 'MultiBatch' or 'Sequential' or 'Gillespie' or 'CRN'. - 'MultiBatch':

    batss_rust.SimulatorMultiBatch does O(sqrt(n)) interaction steps in parallel using batching, and is much faster for large population sizes and relatively small state sets.

    • 'Gillespie':

      uncondtionally uses the Gillespie algorithm. Still uses the multibatch simulator, but instructs it to always use the Gillespie algorithm.

    • 'Sequential':

      batss_rust.SimulatorSequentialArray represents the population as an array of agents, and simulates each interaction step by choosing a pair of agents to update. Defaults to ‘MultiBatch’.

    • 'CRN':

      batss_rust.SimulatorCRNMultiBatch does parallel batching for arbitrary CRNs, and should be faster than Gillespie on large population sizes and small species sets.

  • transition_order (str) –

    Should the rule be interpreted as being symmetric, either 'asymmetric', 'symmetric', or 'symmetric_enforced'. Defaults to ‘symmetric’.

    'asymmetric':

    Ordering of the inputs matters, and all inputs not explicitly given as assumed to be null interactions.

    'symmetric':

    The input pairs are interpreted as unordered. If rule(a, b) returns None, while rule(b, a) has a non-null output, then the output of rule(a, b) is assumed to be the same as rule(b, a). If rule(a, b) and rule(b, a) are each given, there is no effect. Asymmetric interactions can be explicitly included this way.

    'symmetric_enforced':

    The same as symmetric, except that if rule(a, b) and rule(b, a) are non-null and do not give the same set of outputs, a ValueError is raised.

  • seed (int | None) – An optional integer used as the seed for all pseudorandom number generation. Defaults to None.

  • volume (float | None) – If a list of ppsim.crn.Reaction objects is given for a CRN, then the parameter volume can be passed in here. Defaults to None. If None, the volume will be assumed to be the population size n.

  • continuous_time (bool) – Whether continuous time is used. Defaults to False. If a CRN as a list of reactions is passed in, this will be set to True.

  • time_units (str | None) – An optional string given the units that time is in. Defaults to None. This must be a valid string to pass as unit to pandas.to_timedelta.

  • **kwargs

    If rule is a function, any extra function parameters are passed in here, beyond the first two arguments representing the two agents. For example, if rule is defined:

    def rule(sender: int, receiver: int, threshold: int) -> Tuple[int, int]:
        if sender + receiver > threshold:
            return 0, 0
        else:
            return sender, receiver+1
    

    To use a threshold of 20 in each interaction, in the Simulation constructor, use

    sim = Simulation(init_config, rule, threshold=20)
    

  • crn (CRN | None)

simulator: Simulator

An internal Simulator that performs the simulation steps. We probably want to ditch this since we only have one simulator now.

discrete_steps_total: list[int]

Parallel to self.times and self.configs, this is how many total steps were taken up to that time, including steps that simulated passive reactions. This is a cumulative list, i.e., self.discrete_steps_total[i] is the total number of steps taken up to time self.times[i].

discrete_steps_total_last_run: list[int]

Like Simulation.discrete_steps_total, but only since the last call to Simulation.simulator.run().

discrete_steps_no_nulls: list[int]

Parallel to self.times and self.configs, this is how many total steps were taken up to that time, NOT including steps that simulated passive reactions. Cumulative list similarly to self.discrete_steps_total.

discrete_steps_no_nulls_last_run: list[int]

Like Simulation.discrete_steps_total, but only since the last call to Simulation.simulator.run().

simulator_method: str
seed: int | None

The optional integer seed used for rng and inside Rust code.

rng: Generator

A numpy random generator used to sample random variables outside Rust code.

steps_per_time_unit: float

Number of simulated interactions per time unit.

time_units: str | None

The units that time is in.

For options see https://pandas.pydata.org/docs/reference/api/pandas.to_timedelta.html

continuous_time: bool

Whether continuous time is used. The regular discrete time model considers Simulation.steps_per_time_unit steps to be 1 unit of time. The continuous time model is a Poisson process, with expected Simulation.steps_per_time_unit steps per 1 unit of time.

state_list: list[Hashable]

A sorted list of all reachable states.

state_dict: dict[Hashable, int]

Maps states to their integer index to be used in array representations.

column_names: MultiIndex | list[str]

Columns representing all states for pandas dataframe. If the State is a tuple, NamedTuple, or dataclass, this will be a pandas MultiIndex based on the various fields. Otherwise it is list of str(State) for each State.

configs: list[ndarray[tuple[Any, ...], dtype[uint64]]]

A list of all configurations that have been recorded during the simulation, as integer arrays.

times: list[float]

A list of all the corresponding times for configs.

time: float

The current time.

snapshots: list[Snapshot]

A list of ppsim.snapshot.Snapshot objects, which get periodically called during the running of the simulation to give live updates.

rule(a, b)[source]

The rule, as a function of two input states.

Parameters:
Return type:

tuple[Hashable, Hashable] | dict[tuple[Hashable, Hashable], float]

initialize_simulator(config)[source]

Build the data structures necessary to instantiate the simulator class.

Parameters:

config (ndarray[tuple[Any, ...], dtype[uint64]]) – The config array to instantiate the simulator.

Return type:

None

array_from_dict(d)[source]

Convert a configuration dictionary to an array.

Parameters:

d (dict) – A dictionary mapping states to counts.

Return type:

ndarray[tuple[Any, ...], dtype[uint64]]

Returns:

An array giving counts of all states, in the order of self.state_list.

silent()[source]

Check if the configuration is silent.

Return type:

bool

Returns:

True if the configuration is silent, False otherwise.

run(run_until=None, history_interval=1.0, stopping_interval=1.0, timer=True)[source]

Runs the simulation.

Can give a fixed amount of time to run the simulation, or a function that checks the configuration for convergence.

Parameters:
  • run_until (Union[float, Callable[[dict[Hashable, int]], bool], None]) – The stop condition. To run for a fixed amount of time, give a numerical value. To run until a convergence criterion, give a function mapping a configuration (as a dictionary mapping states to counts) to a boolean. The run will stop when the function returns True. Defaults to None. If None, the simulation will run until the configuration is silent (all transitions are null). This only works with the multibatch simulator method, if another simulator method is given, then using None will raise a ValueError.

  • history_interval (Union[float, Callable[[float], float]]) – The length to run the simulator before recording data, in current time units. Defaults to 1. This can either be a float, or a function that takes the current time and and returns a float.

  • stopping_interval (float) – The length to run the simulator before checking for the stop condition.

  • timer (bool) – If True, and there are no other snapshot objects, a default ppsim.snapshot.TimeUpdate snapshot will be created to print updates with the current time. Defaults to True.

Return type:

None

property reactions: str

A string showing all non-null transitions in reaction notation.

Each reaction is separated by newlines, so that print(self.reactions) will display all reactions. Only works with simulator method multibatch, otherwise will raise a ValueError.

property enabled_reactions: str

A string showing all non-null transitions that are currently enabled.

This can only check the current configuration in self.simulator. Each reaction is separated by newlines, so that print(self.enabled_reactions) will display all enabled reactions.

reset(init_config=None)[source]

Reset the simulation.

Parameters:

init_config (dict[Hashable, int] | None) – The configuration to reset to. Defaults to None. If None, will use the old initial configuration.

Return type:

None

set_config(config)[source]

Change the current configuration.

Parameters:

config (dict[Hashable, int] | ndarray) – The configuration to change to. This can be a dictionary, mapping states to counts, or an array giving counts in the order of Simulation.state_list.

Return type:

None

time_to_steps(time)[source]

Convert simulated time into number of simulated interaction steps.

Parameters:

time (float) – The amount of time to convert.

Return type:

int

property config_dict: dict[Hashable, int]

The current configuration, as a dictionary mapping states to counts.

property config_array: ndarray

The current configuration in the simulator, as an array of counts.

The array is given in the same order as self.state_list. The index of state s is self.state_dict[s].

property history: DataFrame

A pandas dataframe containing the history of all recorded configurations.

property null_probability: float

The probability the next interaction is null.

times_in_units(times)[source]

If Simulation.time_units is defined, convert time list to appropriate units.

Parameters:

times (Iterable[float])

Return type:

Iterable[Any]

add_config()[source]

Appends the current simulator configuration and time.

Return type:

None

set_snapshot_time(time)[source]

Updates all snapshots to the nearest recorded configuration to a specified time.

Parameters:

time (float) – The parallel time to update the snapshots to.

Return type:

None

set_snapshot_index(index)[source]

Updates all snapshots to the configuration Simulation.configs [index].

Parameters:

index (int) – The index of the configuration.

Return type:

None

add_snapshot(snap)[source]

Add a new ppsim.snapshot.Snapshot to Simulation.snapshots.

Parameters:

snap (Snapshot) – The ppsim.snapshot.Snapshot object to be added.

Return type:

None

snapshot_slider(var='index')[source]

Returns a slider that updates all ppsim.snapshot.Snapshot objects.

Returns a slider from the ipywidgets library.

Parameters:

var (str) – What variable the slider uses, either 'index' or 'time'.

Return type:

Any

sample_silence_time()[source]

Starts a new trial from the initial distribution and return time until silence.

Return type:

float

sample_future_configuration(time, num_samples=100)[source]

Repeatedly samples the configuration at a fixed future time.

Parameters:
  • time (float) – The amount of time ahead to sample the configuration.

  • num_samples (int) – The number of samples to get.

Return type:

DataFrame

Returns:

A dataframe whose rows are the sampled configuration.

batss.simulation.time_trials(rule, ns, initial_conditions, convergence_condition=None, convergence_check_interval=0.1, num_trials=100, max_wallclock_time=86400, **kwargs)[source]

Gathers data about the convergence time of a rule.

Parameters:
  • rule (Union[Callable[[Hashable, Hashable], tuple[Hashable, Hashable] | dict[tuple[Hashable, Hashable], float]], dict[tuple[Hashable, Hashable], tuple[Hashable, Hashable] | dict[tuple[Hashable, Hashable], float]], Iterable[Reaction]]) – The rule that is used to generate the Simulation object.

  • ns (list[int]) – A list of population sizes n to sample from. This should be in increasing order.

  • initial_conditions (Union[Callable[[int], dict[Hashable, int]], list[dict[Hashable, int]]]) – An initial condition is a dict mapping states to counts. This can either be a list of initial conditions, or a function mapping population size n to an initial condition of n agents.

  • convergence_condition (Optional[Callable[[dict[Hashable, int]], bool]]) – A boolean function that takes a configuration dict (mapping states to counts) as input and returns True if that configuration has converged. Defaults to None. If None, the simulation will run until silent (all transitions are null), and the data will be for silence time.

  • convergence_check_interval (float) – How often (in parallel time) the simulation will run between convergence checks. Defaults to 0.1. Smaller values give better resolution, but spend more time checking for convergence.

  • num_trials (int) – The maximum number of trials that will be done for each population size n, if there is sufficient time. Defaults to 100. If you want to ensure that you get the full num_trials samples for each value n, use a large value for time_bound.

  • max_wallclock_time (float) – A bound (in seconds) for how long this function will run. Each value n is given a time budget based on the time remaining, and will stop before doing num_trials runs when this time budget runs out. Defaults to 60 * 60 * 24 (one day).

  • **kwargs – Other keyword arguments to pass into Simulation.

Return type:

DataFrame

Returns:

A pandas dataframe giving the data from each trial, with two columns 'n'' and 'time'. A good way to visualize this dataframe is using the seaborn library, calling sns.lineplot(x='n', y='time', data=df).

Module for expressing population protocols using CRN notation. Ideas and much code taken from https://github.com/enricozb/python-crn.

The general syntax is

a, b, u = species('A B U')
approx_majority = [
    a + b >> 2 * u,
    a + u >> 2 * a,
    b + u >> 2 * b,
]
n = 10 ** 5
init_config = {a: 0.51 * n, b: 0.49 * n}
sim = Simulation(init_config=init_config, rule=approx_majority)

In other words, a list of reactions is treated by the batss library just like the other ways of specifying population protocol transitions (the rule parameter in the constructor for batss.simulation.Simulation, which also accepts a dict or a Python function).

More examples given in https://github.com/UC-Davis-molecular-computing/batss/tree/main/examples

batss.crn.species(sp)[source]

Create a list of Specie (Single species Expression’s), or a single one.

Parameters:

sp (Union[str, Iterable[str]]) – An string or Iterable of strings representing the names of the species being created. If a single string, species names are interpreted as space-separated.

Return type:

tuple[Specie, ...]

Examples:

w, x, y, z = species('W X Y Z')
rxn = x + y >> z + w
w, x, y, z = species(['W', 'X', 'Y', 'Z'])
rxn = x + y >> z + w
batss.crn.replace_reversible_rxns(rxns)[source]
Parameters:

rxns (Iterable[Reaction]) – list of Reaction’s

Return type:

list[Reaction]

Returns: list of Reaction’s, where every reversible reaction in rxns has been replaced by

two irreversible reactions, and all others have been left as they are

batss.crn.reactions_to_dict(reactions, n, volume)[source]

Returns dict representation of reactions, transforming unimolecular reactions to bimolecular, and converting rates to probabilities, also returning the max rate so the ppsim.simulation.Simulation knows how to scale time.

Parameters:
  • reactions (Iterable[Reaction]) – list of Reaction’s

  • n (int) – the population size, necessary for rate conversion

  • volume (float) – parameter as defined in Gillespie algorithm

Return type:

tuple[dict[tuple[Specie, Specie], tuple[Specie, Specie] | dict[tuple[Specie, Specie], float]], float]

Returns:

(transitions_dict, max_rate), where transitions_dict is the dict representation of the transitions, and max_rate is the maximum rate for any pair of reactants, i.e., if we have reactions (a + b >> c + d).k(2) and (a + b >> x + y).k(3), then the ordered pair (a,b) has rate 2+3 = 5

batss.crn.convert_unimolecular_to_bimolecular_and_flip_reactant_order(reactions, n, volume)[source]

Process all reactions before being added to the dictionary.

bimolecular reactions have their rates multiplied by the corrective factor (n-1) / (2 * volume). Bimolecular reactions with two different reactants are added twice, with their reactants in both orders.

Parameters:
Return type:

list[Reaction]

class batss.crn.Specie(name, is_special_specie=False)[source]

Bases: object

Parameters:
  • name (str)

  • is_special_specie (bool)

name: str
is_special_specie: bool = False
class batss.crn.Expression(species)[source]

Bases: object

Class used for very basic symbolic manipulation of left/right hand side of stoichiometric equations. Not very user friendly; users should just use the species functions and manipulate those to get their reactions.

Parameters:

species (list[Specie])

species: list[Specie]

ordered list of species in expression, e.g, A+A+B is [A,A,B]

get_species()[source]

Returns the set of species in this expression, not their coefficients.

Return type:

set[Specie]

batss.crn.concentration_to_count(concentration, volume)[source]
Parameters:
  • concentration (float) – units of M (molar) = moles / liter

  • volume (float) – units of liter

Return type:

int

Returns:

count of molecule with concentration in volume

class batss.crn.RateConstantUnits(*values)[source]

Bases: Enum

stochastic = 'stochastic'

Units of L/s. Multiple by Avogadro’s number to convert to mass-action units.

mass_action = 'mass_action'

Units of /M/s. Divide by Avogadro’s number to convert to stochastic units.

class batss.crn.Reaction(reactants, products, k=1, r=1, rate_constant_units=RateConstantUnits.stochastic, rate_constant_reverse_units=RateConstantUnits.stochastic, reversible=False)[source]

Bases: object

Representation of a stoichiometric reaction using a pair of Expressions, one for the reactants and one for the products.

Parameters:
  • reactants (Specie | Expression)

  • products (Specie | Expression)

  • k (float)

  • r (float)

  • rate_constant_units (RateConstantUnits)

  • rate_constant_reverse_units (RateConstantUnits)

  • reversible (bool)

  • reactants – left side of species in the reaction

  • products – right side of species in the reaction

  • k – Rate constant of forward reaction

  • r – Rate constant of reverse reaction (only used if Reaction.reversible is true

  • rate_constant_units – Units of forward rate constant

  • rate_constant_reverse_units – Units of reverse rate constant

  • reversible – Whether reaction is reversible

reactants: Expression

The left side of species in the reaction.

products: Expression

The right side of species in the reaction.

rate_constant: float = 1

Rate constant of forward reaction.

rate_constant_reverse: float = 1

Rate constant of reverse reaction (only used if Reaction.reversible is true).

rate_constant_units: RateConstantUnits = 'stochastic'

Units of forward rate constant.

rate_constant_reverse_units: RateConstantUnits = 'stochastic'

Units of reverse rate constant.

reversible: bool = False

Whether reaction is reversible, i.e. products → reactants is a reaction also.

is_unimolecular()[source]

Returns: true if there is one reactant

Return type:

bool

is_bimolecular()[source]

Returns: true if there are two reactants

Return type:

bool

symmetric()[source]

Returns: true if there are two reactants that are the same species

Return type:

bool

symmetric_products()[source]

Returns: true if there are two products that are the same species

Return type:

bool

num_reactants()[source]

Returns: number of reactants

Return type:

int

num_products()[source]

Returns: number of products

Return type:

int

is_conservative()[source]

Returns: true if number of reactants equals number of products

Return type:

bool

reactant_if_unimolecular()[source]

Returns: unique reactant if there is only one Raises: ValueError if there are multiple reactants

Return type:

Specie

product_if_unique()[source]

Returns: unique product if there is only one Raises: ValueError if there are multiple products

Return type:

Specie

reactants_if_bimolecular()[source]

Returns: pair of reactants if there are exactly two Raises: ValueError if there are not exactly two reactants

Return type:

tuple[Specie, Specie]

reactant_names_if_bimolecular()[source]

Returns: pair of reactant names if there are exactly two Raises: ValueError if there are not exactly two reactants

Return type:

tuple[str, str]

products_if_exactly_two()[source]

Returns: pair of products if there are exactly two Raises: ValueError if there are not exactly two products

Return type:

tuple[Specie, Specie]

product_names_if_exactly_two()[source]

Returns: pair of product names if there are exactly two Raises: ValueError if there are not exactly two products

Return type:

tuple[str, str]

generativity()[source]

Returns: number of products minus number of reactants

Return type:

int

property rate_constant_stochastic: float

forward rate constant in stochastic units (converts from mass-action if necessary)

Type:

Returns

property rate_constant_reverse_stochastic: float

reverse rate constant in stochastic units (converts from mass-action if necessary)

Type:

Returns

k(coeff, units=RateConstantUnits.stochastic)[source]

Changes the reaction coefficient to coeff and returns self.

This is useful for including the rate constant during the construction of a reaction. For example

x, y, z = species("X Y Z")
rxns = [
    (x + y >> z).k(2.5),
    (z >> x).k(1.5),
    (z >> y).k(0.5)),
]
Parameters:
  • coeff (float) – float The new reaction coefficient

  • units (RateConstantUnits) – float units of rate constant (default stochastic)

Return type:

Reaction

r(coeff, units=RateConstantUnits.stochastic)[source]

Changes the reverse reaction reaction rate constant to coeff and returns self.

This is useful for including the rate constant during the construction of a reaction. For example

x, y, z = species("X Y Z")
rxns = [
    (x + y >> z).k(2.5),
    (z >> x).k(1.5),
    (z >> y).k(0.5)),
]
Parameters:
  • coeff (float) – float The new reverse reaction rate constant

  • units (RateConstantUnits) – float units of rate constant (default stochastic)

Return type:

Reaction

get_species()[source]
Return type:

tuple[Specie]

Returns:

the species present in the products and reactants.

class batss.crn.CRN(reactions, species)[source]

Bases: object

Representation of a CRN, a set of reactions.

Parameters:
reactions: list[Reaction]

The reactions comprising the CRN.

species: list[Specie]

The species in the CRN. Ordered by index that represents them during algorithm execution.

order()[source]

Returns: the order of the CRN, the greatest number of reactants in any reaction

Return type:

int

generativity()[source]

Returns: the generativity of the CRN, the greatest number of products minus reactants in any reaction

Return type:

int

exportable_reactions()[source]

Returns: the reactions of the CRN formatted in a way that can be passed through pyo3.

Return type:

list[tuple[list[int], list[int], float]]

batss.crn.extra_species()[source]
batss.crn.convert_to_uniform(crn, volume)[source]

Convert a CRN to an equivalent uniform CRN. The new CRN will have two new species, K (catalyst) and W (waste). This function does NOT adjust rate constants for the count of K; this is done later.

Parameters:
  • crn (CRN) – a CRN to be converted. Should not be an output of this function.

  • volume (float) – volume of the CRN.

Return type:

CRN

batss.crn.catalyst_specie()[source]
Return type:

Specie

batss.crn.waste_specie()[source]
Return type:

Specie

batss.crn.species_in_rxns(rxns)[source]
Parameters:

rxns (Iterable[Reaction]) – iterable of Reaction’s

Return type:

list[Specie]

Returns: list of species (without repetitions) in Reaction’s in rxns

batss.crn.gpac_format(rxns, init_config)[source]

Create a gpac CRN in the form of equivalent initial configuration and list of gpac Reaction objects.

These can be simulated using the Gillespie algorithm (with the rebop package as a backend) in gpac using gpac.rebop_crn_counts, and plotted using gpac.plot_crn_counts.

Parameters:
Return type:

tuple[list[Reaction], dict[Specie, int]]

Returns:

A tuple of (reactions, config) essentially equivalent to the ppsim init_config and rxns but using the gpac package’s versions of those data structures.

batss.crn.gillespy2_format(init_config, rxns, volume=1.0)[source]

Create a gillespy2 Model object from a CRN description.

Parameters:
Return type:

Any

Returns:

An equivalent gillespy2 Model object

batss.crn.stochkit_format(rxns, init_config, volume=1.0, name='CRN')[source]
Parameters:
  • rxns (Iterable[Reaction]) – reactions to translate to StochKit format

  • init_config (dict[Specie, int]) – dict mapping each Specie to its initial count

  • volume (float) – volume in liters

  • name (str) – name of the CRN

Return type:

str

Returns:

string describing CRN in StochKit XML format

batss.crn.write_stochkit_file(filename, rxns, init_config, volume=1.0, name='CRN')[source]

Write stochkit file :type filename: str :param filename: name of file to write :type rxns: Iterable[Reaction] :param rxns: reactions to translate to StochKit format :type init_config: dict[Specie, int] :param init_config: dict mapping each Specie to its initial count :type volume: float :param volume: volume in liters :type name: str :param name: name of the CRN

Return type:

None

Parameters:

A module for Snapshot objects used to visualize the protocol during or after the simulation has run.

Snapshot is a base class for snapshot objects that get are updated by batss.simulation.Simulation.

Plotter is a subclass of Snapshot that creates a matplotlib figure and axis. It also gives the option for a state_map function which maps states to the categories which will show up in the plot.

StatePlotter is a subclass of Plotter that creates a barplot of the counts in categories.

HistoryPlotter is a subclass of Plotter that creates a lineplot of the counts in categories over time.

class batss.snapshot.Snapshot(update_time=0.1)[source]

Bases: ABC

Abstract base class for snapshot objects.

Parameters:

update_time (float)

Init constructor for the base class.

Parameters can be passed in here, and any attributes that can be defined without the parent ppsim.simulation.Simulation object can be instantiated here, such as Snapshot.update_time.

Parameters:

update_time (float)

simulation: Optional[Simulation]

The ppsim.simulation.Simulation object that initialized and will update the Snapshot. This attribute gets set when the ppsim.simulation.Simulation object calls ppsim.simulation.Simulation.add_snapshot().

update_time: float

How many seconds will elapse between calls to update while in the ppsim.simulation.Simulation.run() method of ppsim.simulation.Simulation.

time: float

The time at the current snapshot. Changes when Snapshot.update() is called.

config: ndarray | None

The configuration array at the current snapshot. Changes when Snapshot.update() is called.

next_snapshot_time: float

The time at which the next snapshot will be taken.

initialize()[source]

Method which is called once during ppsim.simulation.Simulation.add_snapshot(). Should only be called after ppsim.simulation.Simulation.add_snapshot() is called.

Any initialization that requires accessing the data in Snapshot.simulation should go here.

Return type:

None

update(index=None)[source]

Method which is called while Snapshot.simulation is running.

Parameters:

index (int | None) – An optional integer index. If present, the snapshot will use the data from configuration ppsim.simulation.Simulation.configs [index] and time ppsim.simulation.Simulation.times [index]. Otherwise, the snapshot will use the current configuration ppsim.simulation.Simulation.config_array() and current time.

Return type:

None

class batss.snapshot.TimeUpdate(time_bound=None, update_time=0.2)[source]

Bases: Snapshot

Simple Snapshot that prints the current time in the ppsim.simulation.Simulation.

When calling ppsim.simulation.Simulation.run, if ppsim.simulation.Simulation.snapshots is empty, then this object will get added to provide a basic progress update.

Parameters:

Init constructor for the base class.

Parameters can be passed in here, and any attributes that can be defined without the parent ppsim.simulation.Simulation object can be instantiated here, such as Snapshot.update_time.

Parameters:
pbar: tqdm

The progress bar object that gets created.

start_time: float

The starting time of the simulation.

initialize()[source]

Method which is called once during ppsim.simulation.Simulation.add_snapshot(). Should only be called after ppsim.simulation.Simulation.add_snapshot() is called.

Any initialization that requires accessing the data in Snapshot.simulation should go here.

Return type:

None

update(index=None)[source]

Method which is called while Snapshot.simulation is running.

Parameters:

index (int | None) – An optional integer index. If present, the snapshot will use the data from configuration ppsim.simulation.Simulation.configs [index] and time ppsim.simulation.Simulation.times [index]. Otherwise, the snapshot will use the current configuration ppsim.simulation.Simulation.config_array() and current time.

Return type:

None

class batss.snapshot.Plotter(state_map=None, update_time=0.5, yscale='linear')[source]

Bases: Snapshot

Base class for a Snapshot which will make a plot.

Gives the option to map states to categories, for an easy way to visualize relevant subsets of the states rather than the whole state set. These require an interactive matplotlib backend to work.

Parameters:

Initializes the Plotter.

Parameters:
  • state_map (Optional[Callable[[Hashable], Any]]) – An optional function mapping states to categories. If None, then the state itself will be used as the category.

  • update_time (float) – How many seconds will elapse between calls to update while ppsim.simulation.Simulation.run method.

  • yscale (str) – The scale used for the yaxis, passed into ax.set_yscale. Defaults to ‘linear’.

fig: Optional[Figure]

The matplotlib figure that is created.

ax: Optional[Axes]

The matplotlib axis object that is created. Modifying properties of this object is the most direct way to modify the plot.

categories: list[Any]

A list which holds the set {state_map(state)} for all states in Snapshot.simulation.state_list.

state_map: Callable[[Hashable], Any]

A function mapping states to categories, which acts as a filter to view a subset of the states or just one field of the states. This function should take a state as input and return a category. If not specified in constructor, then the state itself will be used as the category.

yscale: str

The scale used for the yaxis, passed into ax.set_yscale.

initialize()[source]

Initializes the plotter by creating a fig and ax.

Return type:

None

class batss.snapshot.StatePlotter(state_map=None, update_time=0.5, yscale='linear')[source]

Bases: Plotter

Plotter which produces a barplot of counts.

Parameters:

Initializes the Plotter.

Parameters:
  • state_map (Optional[Callable[[Hashable], Any]]) – An optional function mapping states to categories. If None, then the state itself will be used as the category.

  • update_time (float) – How many seconds will elapse between calls to update while ppsim.simulation.Simulation.run method.

  • yscale (str) – The scale used for the yaxis, passed into ax.set_yscale. Defaults to ‘linear’.

initialize()[source]

Initializes the barplot.

Return type:

None

If Plotter.state_map gets changed, call Plotter.initialize() to update the barplot to

show the new set Plotter.categories.

update(index=None)[source]

Update the heights of all bars in the plot.

Parameters:

index (int | None)

Return type:

None

class batss.snapshot.HistoryPlotter(state_map=None, update_time=0.5, yscale='linear')[source]

Bases: Plotter

Plotter which produces a lineplot of counts over time.

Parameters:

Initializes the Plotter.

Parameters:
  • state_map (Optional[Callable[[Hashable], Any]]) – An optional function mapping states to categories. If None, then the state itself will be used as the category.

  • update_time (float) – How many seconds will elapse between calls to update while ppsim.simulation.Simulation.run method.

  • yscale (str) – The scale used for the yaxis, passed into ax.set_yscale. Defaults to ‘linear’.

update(index=None)[source]

Make a new history plot.

Parameters:

index (int | None)

Return type:

None