“... even if the previous millisecond is closer to us than the birth of the universe, it is equally out of reach.” ― Jean-Christophe Valtat, Luminous Chaos
This library implements a generalized version of the Gillespie
Algorithm, a stochastic
approach to numerically solving discrete systems. Each iteration, the algorithm
will calculate the propensities for each reaction given a rate and the counts
of the reactants present in the current state of the system, then selects one
reaction to occur and the interval of time between the previous reaction and
the current reaction. Iterating this produces a trajectory (or history
) of
the state vector over the course of the simulation.
Add the following to your requirements.txt
, or run
pip install stochastic-arrow
to install it from PyPI:
stochastic-arrow
NOTE: If upgrading from a version older than 1.0.0, check if the arrow
datetime package is installed. If so, uninstall arrow
before upgrading stochastic-arrow
, then reinstall arrow
.
> pip show arrow
> pip uninstall arrow
> pip install stochastic-arrow
> pip install arrow
The stochastic_arrow
library presents a single class as an interface,
StochasticSystem
, which operates on a set of reactions (encoded as a numpy
matrix of stoichiometrix coefficients) and associated reaction rates:
from stochastic_arrow import StochasticSystem
import numpy as np
# Each row is a reaction and each column is a molecular species (or other
# entity). The first reaction here means that the first and second elements
# combine to create the third, while the fourth is unaffected.
stoichiometric_matrix = np.array([
[1, 1, -1, 0],
[-2, 0, 0, 1],
[-1, -1, 1, 0]], np.int64)
# Once we have a matrix of reactions, we can
# construct the system.
system = StochasticSystem(stoichiometric_matrix)
Now that the system has been instantiated, we can invoke it with any initial state vector and set of reaction rates and then run it for a given time interval:
# This gives the initial state of the system (counts of each molecular species,
# for instance).
import numpy as np
state = np.array([1000, 1000, 0, 0])
# We also specify how long we want the simulation to run. Here we set it to one
# second.
duration = 1
# Each reaction has an associated rate for how probable that reaction is.
rates = np.array([3.0, 1.0, 1.0])
Once we have an initial state and rates, we can run the simulation for the
given duration. evolve
returns a dictionary with five keys:
- steps - the number of steps the simulation took
- time - at what time point each event took place
- events - the events that occurred
- occurrences - the number of times each event occurred (derived directly from
events
) - outcome - the final state of the system
result = system.evolve(duration, state, rates)
If you are interested in the history of states for plotting or otherwise, these can be
derived from the list of events and the stoichiometric matrix, along with the inital
state. reenact_events
will do this for you:
from stochastic_arrow import reenact_events
history = reenact_events(stoichiometric_matrix, result['events'], state)
stochastic_arrow
uses pytest. To test it:
> make clean compile
> pytest
NOTE: make compile
without an explicit clean
might not fully build the extension.
There are more command line features in test_arrow:
> python -m stochastic_arrow.test.test_arrow --complexation
> python -m stochastic_arrow.test.test_arrow --plot
> python -m stochastic_arrow.test.test_arrow --obsidian
> python -m stochastic_arrow.test.test_arrow --memory
> python -m stochastic_arrow.test.test_arrow --time
More examples:
> python -m stochastic_arrow.test.test_hang
> pytest -m stochastic_arrow/test/test_arrow.py
> pytest -k flagella
- Rename module to
stochastic_arrow
to avoid name conflict (Issue #51). All users must update their import statements to usestochastic_arrow
instead ofarrow
.
- Update to Cython 0.29.34. (Cython 3.0.0 is now in beta.)
- Update to Cython 3.0.0a11 for compatibility with Python 3.11.
Add
arrow.pxd
to work around a Cython 3.0.0 bug. - Stop using deprecated
numpy.distutils
to avoid warnings and prepare for its removal in Python 3.12. - Make
test_arrow.py --plot
compatible with Python 3. - Fix
PytestReturnNotNoneWarning
warnings from pytest 7.2.0.
- Add the arrow_hang unit test which catches a nasty edge-case (Issue #48), fix the bug, and make the code more robust to some other potential bugs.
- Can pickle StochasticSystem instances.
- Introduced backwards-incompatible API change for supplying rates at
evolve()
time rather than__init__()
forStochasticSystem
.