chain_simulator.simulation module#

Core simulation functionality.

This module provides the core functionality for simulating digital twins using Markov Chains. There are both high- and low-level functions available. The high-level functions automatically makes use of optimal low-level functions, based on input parameters. The low-level functions are also available in case the high-level functions aren’t of any use.

chain_simulator.simulation.state_vector_processor(state_vector, transition_matrix, steps, interval=None)#

Simulate a Markov chain and return intermediary/final state vector(s).

Dynamically simulate a Markov chain on either a Central Processing Unit (CPU) or Graphics Processing Unit (GPU). The state_vector is multiplied with the transition_matrix, returning an intermediate/final state vector of the n-th or step-th step in time. By default, only a final state vector is returned. Intermediate state vectors are obtained by setting interval, which will represent every n-th intermediate state vector.

Parameters:

state_vectorany 1D array: A 1D array with an initial state probability distribution, i.e. an initial state vector.
transition_matrixany 2D array: A 2D array with state change probabilities, i.e. a transition matrix.
stepsint: How many steps in time transition_matrix must progress.
intervalint, optional: Which n-th or interval-th intermediate state vector must be returned, none by default.

Yields:

tuple of array and int: An intermediate/final state vector of the current step in time and the current step in time.

Raises:

TypeError: If the transition matrix type is incompatible with a GPU.

Parameters:

state_vector (STATE_VECTOR) –
transition_matrix (TRANSITION_MATRIX) –
steps (int) –
interval (Optional[int]) –

Return type:

Iterator[Tuple[NDArray[Any], int]]

See also

chain_simulator: Progress a Markov chain forward in time.

Notes

There are three distinct implementations. The first implementation is GPU-based, implemented with CuPy. CuPy is however an optional dependency. If the library is not installed, or the state vector / transition matrix format is incompatible with a GPU, the function falls back to a SciPy implementation. If the transition matrix is a NumPy array, the function falls back to the final implementation, which is implemented using NumPy.

Vector-matrix multiplication is used with the dot-method. CuPy, SciPy and NumPy all have their own optimised version of this method, so it is important to choose the right implementation. Choosing a wrong implementation results in errors (types are incompatible) and lower performance (a dot product may be faster on a GPU than a CPU).

Examples

Simulate a Markov chain for 3 days with a NumPy transition matrix:

>>> import numpy as np
>>> from chain_simulator.simulation import state_vector_processor
>>> initial_state_vector = np.array([1, 0, 0])
>>> transition_matrix = np.array(
...     [[0.0, 1.0, 0.0], [0.0, 0.5, 0.5], [0.0, 0.0, 1.0]]
... )
>>> simulator = state_vector_processor(
...     initial_state_vector, transition_matrix, 3
... )
>>> next(simulator)
(array([0.  , 0.25, 0.75]), 3)

Simulate a Markov chain for 2 days with a SciPy transition matrix and all intermediary results:

>>> from scipy import sparse
>>> csc_transition_matrix = sparse.csc_array(transition_matrix)
>>> simulator = state_vector_processor(
...     initial_state_vector, csc_transition_matrix, 2, interval=1
... )
>>> next(simulator)
(array([0., 1., 0.]), 1)
>>> next(simulator)
(array([0. , 0.5, 0.5]), 2)

chain_simulator.simulation.vector_processor_cupy(state_vector, transition_matrix, steps, interval=None)#

Process state vectors using CuPy.

Process a state vector using a CuPy implementation. CuPy adds GPU support for processing transition_matrix in both dense and sparse formats. This function multiplies state_vector with transition_matrix to get intermediate/final state vectors.

Parameters:

state_vector1D cupy array: A 1D array with an initial state probability distribution, i.e. an initial state vector.
transition_matrix2D cupy dense/sparse array: A 2D array with state change probabilities, i.e. a transition matrix.
stepsint: How many steps in time transition_matrix must progress.
intervalint, optional: Which n-th or interval-th intermediate state vector must be returned, none by default.

Yields:

tuple of array and int: An intermediate/final state vector of the current step in time and the current step in time.

Raises:

TypeError: If state_vector is not a CuPy array or transition_matrix is not a CuPy dense/sparse array.

Parameters:

state_vector (_cupy.ndarray) –
transition_matrix (CUPY_MATRIX) –
steps (int) –
interval (Optional[int]) –

Return type:

Iterator[Tuple[_cupy.ndarray, int]]

See also

chain_simulator: Progress a Markov chain forward in time.
vector_processor_numpy: Process state vectors using NumPy.
vector_processor_scipy: Process state vectors using SciPy.
state_vector_processor

Notes

When a CPU implementation is not fast enough ( vector_processor_numpy() or vector_processor_scipy()) and a GPU is available, it is best to process a transition matrix on a GPU. GPU support is enabled by installing the optional dependency CuPy and allows processing of both dense and sparse matrices.

It is important that both the state vector and transition matrix are already on the GPU, i.e. cupy.ndarray or cupyx.scipy.sparse.csc_matrix. If one of them is not on the GPU, processing will fail.

This function is slightly lower-level than state_vector_processor() as it does not perform any type conversions. The function chain_simulator() is used to progress the transition matrix forward in time.

Examples

>>> import cupy
>>> initial_state_vector = cupy.array([1, 0, 0])
>>> transition_matrix = cupy.array(
...     [[0.0, 1.0, 0.0], [0.0, 0.5, 0.5], [0.0, 0.0, 1.0]]
... )
>>> simulator = vector_processor_cupy(
...     initial_state_vector, transition_matrix, 3
... )
>>> next(simulator)
(array([[0, 1 / 4, 3 / 4], [0, 1 / 8, 7 / 8], [0, 0, 1]]), 3)

>>> import cupyx.scipy.sparse
>>> csr_transition_matrix = sparse.csr_array(initial_state_vector)
>>> simulator = vector_processor_cupy(
...     initial_state_vector, csr_transition_matrix, 2, interval=1
... )
>>> next(simulator)
(array([[0, 1, 0], [0, 1 / 2, 1 / 2], [0, 0, 1]]), 1)
>>> next(simulator)
(array([[0, 1 / 2, 1 / 2], [0, 1 / 4, 3 / 4], [0, 0, 1]]), 2)

chain_simulator.simulation.vector_processor_numpy(state_vector, transition_matrix, steps, interval=None)#

Process state vectors using NumPy.

Process a state vector using a NumPy implementation. This function multiplies state_vector with transition_matrix to get intermediate/ final state vectors.

Parameters:

state_vector1D numpy array: A 1D array with an initial state probability distribution, i.e. an initial state vector.
transition_matrix2D numpy array: A 2D array with state change probabilities, i.e. a transition matrix.
stepsint: How many steps in time transition_matrix must progress.
intervalint, optional: Which n-th or interval-th intermediate state vector must be returned, none by default.

Yields:

tuple of array and int: An intermediate/final state vector of the current step in time and the current step in time.

Raises:

TypeError: If both state vector and transition matrix are not instances of numpy.ndarray.

Parameters:

state_vector (NDArray[Any]) –
transition_matrix (NDArray[Any]) –
steps (int) –
interval (Optional[int]) –

Return type:

Iterator[Tuple[NDArray[Any], int]]

See also

chain_simulator: Progress a Markov chain forward in time.
vector_processor_scipy: Process state vectors using SciPy.
vector_processor_cupy: Process state vectors using CuPy.
state_vector_processor

Notes

When a transition matrix is not all that sparse (density > 33%) or it is small enough to fit in memory, it is best to process them using a NumPy implementation. This way there is no overhead of sparse format conversions or sparse format storage.

Examples

>>> import numpy as np
>>> initial_state_vector = np.array([1, 0, 0])
>>> transition_matrix = np.array(
...     [[0.0, 1.0, 0.0], [0.0, 0.5, 0.5], [0.0, 0.0, 1.0]]
... )
>>> simulator = vector_processor_numpy(
...     initial_state_vector, transition_matrix, 3
... )
>>> next(simulator)
(array([0.  , 0.25, 0.75]), 3)

>>> simulator = vector_processor_numpy(
...     initial_state_vector, transition_matrix, 2, interval=1
... )
>>> next(simulator)
(array([0., 1., 0.]), 1)
>>> next(simulator)
(array([0. , 0.5, 0.5]), 2)

chain_simulator.simulation.vector_processor_scipy(state_vector, transition_matrix, steps, interval=None)#

Process state vectors using SciPy.

Process a state vector using a SciPy implementation. SciPy adds support for processing transition_matrix in sparse formats. This function multiplies state_vector with transition_matrix to get intermediate/final state vectors.

Parameters:

state_vector1D numpy array: A 1D array with an initial state probability distribution, i.e. an initial state vector.
transition_matrix2D scipy sparse array/matrix: A 2D array with state change probabilities, i.e. a transition matrix.
stepsint: How many steps in time transition_matrix must progress.
intervalint, optional: Which n-th or interval-th intermediate state vector must be returned, none by default.

Yields:

tuple of array and int: An intermediate/final state vector of the current step in time and the current step in time.

Raises:

TypeError: If state_vector is not of type numpy.ndarray or transition_matrix is not a SciPy CSC/CSR array/matrix.

Parameters:

state_vector (NDArray[Any]) –
transition_matrix (SCIPY_SPARSE_MATRIX) –
steps (int) –
interval (Optional[int]) –

Return type:

Iterator[Tuple[NDArray[Any], int]]

See also

chain_simulator: Progress a Markov chain forward in time.
vector_processor_numpy: Process state vectors using NumPy.
vector_processor_cupy: Process state vectors using CuPy.
state_vector_processor

Notes

When a transition matrix no longer fits in memory as a dense format, sparse formats from scipy. sparse are available. These sparse formats provide their own vector-matrix multiplication functionality. If this is not used, a sparse transition matrix is converted into a regular/dense NumPy matrix and fits most likely no longer in memory.

Not all sparse formats can be used for arithmetic operations. Compressed Sparse Column (CSC) and Compressed Sparse Row (CSR) formats allow for efficient arithmetic operations, while COOrdinate (COO) format does not support any arithmetic operations.

Examples

>>> import numpy as np
>>> import scipy
>>> initial_state_vector = np.array([1, 0, 0])
>>> transition_matrix = scipy.sparse.csc_array(
...     [[0.0, 1.0, 0.0], [0.0, 0.5, 0.5], [0.0, 0.0, 1.0]]
... )
>>> simulator = vector_processor_scipy(
...     initial_state_vector, transition_matrix, 3
... )
>>> next(simulator)
(array([0.  , 0.25, 0.75]), 3)

>>> simulator = vector_processor_scipy(
...     initial_state_vector, transition_matrix, 2, interval=1
... )
>>> next(simulator)
(array([0., 1., 0.]), 1)
>>> next(simulator)
(array([0. , 0.5, 0.5]), 2)