stochastic
Module/class that carries out different type of simulation on an ode formulation
-
class
pygom.model.simulate.
SimulateOde
(state=None, param=None, derived_param=None, transition=None, birth_death=None, ode=None) This builds on top of
DeterministicOde
which we simulate the outcome instead of solving it deterministicallyParameters: - state: list
A list of states (string)
- param: list
A list of the parameters (string)
- derived_param: list
A list of the derived parameters (tuple of (string,string))
- transition: list
A list of transition (
Transition
)- birth_death: list
A list of birth or death process (
Transition
)- ode: list
A list of ode (
Transition
)
-
birth_death_rate
(state, t) Evaluate the birth death rates given state and time
Parameters: - state: array like
The current numerical value for the states which can be
np.ndarray
orlist
- t: double
The current time
Returns: numpy.ndarray
an array of size (M,M) where M is the number of birth and death actions
-
cle
(x0, t0, t1, output_time=False) Stochastic simulation using the CLE approximation starting from time t0 to t1 with the starting state values of x0. The CLE approximation is performed using a simple Euler-Maruyama method with step size h. We assume that the input parameter transition_func provides \(f(x,t)\) while the CLE is defined as \(dx = x + V*h*f(x,t) + \sqrt(f(x,t))*Z*\sqrt(h)\) with \(Z\) being standard normal random variables.
Parameters: - x: array like
state vector
- t0: double
start time
- t1: double
final time
-
eval_birth_death_rate
(parameters=None, time=None, state=None) Evaluate the birth and death rates given parameters, state and time.
Parameters: - parameters: list
see
setParameters()
- time: double
The current time
- state: array list
The current numerical value for the states which can be
numpy.ndarray
orlist
Returns: numpy.matrix
ormpmath.matrix
Matrix of dimension [number of birth and death rates x 1]
-
eval_transition_matrix
(parameters=None, time=None, state=None) Evaluate the transition matrix given parameters, state and time. Note that the output is not in sparse format
Parameters: - parameters: list
see
setParameters()
- time: double
The current time
- state: array list
The current numerical value for the states which can
numpy.ndarray
orlist
Returns: numpy.matrix
ormpmath.matrix
Matrix of dimension [number of state x number of state]
-
eval_transition_mean
(parameters=None, time=None, state=None) Evaluate the transition mean given parameters, state and time.
Parameters: - parameters: list
see
setParameters()
- time: double
The current time
- state: array list
The current numerical value for the states which can be
numpy.ndarray
orlist
Returns: numpy.matrix
ormpmath.matrix
Matrix of dimension [number of state x number of state]
-
eval_transition_var
(parameters=None, time=None, state=None) Evaluate the transition variance given parameters, state and time.
Parameters: - parameters: list
see
setParameters()
- time: double
The current time
- state: array list
The current numerical value for the states which can be
numpy.ndarray
orlist
Returns: numpy.matrix
ormpmath.matrix
Matrix of dimension [number of state x number of state]
-
eval_transition_vector
(parameters=None, time=None, state=None) Evaluate the transition vector given parameters, state and time. Note that the output is not in sparse format
Parameters: - parameters: list
see
setParameters()
- time: double
The current time
- state: array list
The current numerical value for the states which can
numpy.ndarray
orlist
Returns: numpy.matrix
ormpmath.matrix
vector of dimension [total number of transitions]
-
exact
(x0, t0, t1, output_time=False) Stochastic simulation using an exact method starting from time t0 to t1 with the starting state values of x0
Parameters: - x: array like
state vector
- t0: double
start time
- t1: double
final time
-
get_bd_from_ode
(A=None) Returns a list of:class:Transition from this object by unrolling the odes. All the elements are of TransitionType.B or TransitionType.D
-
get_birth_death_rate
() Find the algebraic equations of birth and death processes
Returns: sympy.matrices.matrices
birth death process in matrix form
-
get_transition_matrix
() Returns the transition matrix in algebraic form.
Returns: sympy.matrices.matrices
A matrix of dimension [number of state x number of state]
-
get_transition_vector
() Returns the set of transitions in a single vector, transitions between state to state first then the birth and death process
Returns: sympy.matrices.matrices
A matrix of dimension [total number of transitions x 1]
-
get_transitions_from_ode
() Returns a list of
Transition
from this object by unrolling the odes. All the elements are of TransitionType.T
-
get_unrolled_obj
() Returns a
SimulateOde
with the same state and parameters as the current object but with the equations defined by a set of transitions and birth death process instead of say, odes
-
hybrid
(x0, t0, t1, output_time=False) Stochastic simulation using an hybrid method that uses either the first reaction method or the \(\tau\)-leap depending on the size of the states and transition rates. Starting from time t0 to t1 with the starting state values of x0.
Parameters: - x: array like
state vector
- t0: double
start time
- t1: double
final time
-
plot
(sim_X=None, sim_T=None) Plot the results of a simulation
Takes the output of a function like simulate_jump
Parameters: - sim_X: list
of length iteration each with (len(t),len(state)) if t is a vector, else it outputs unequal shape that was record of all the jumps
- sim_T: list or :class:`numpy.ndarray`
if t is a single value, it outputs unequal shape that was record of all the jumps. if t is a vector, it outputs t so that it is a
numpy.ndarray
instead
Notes
If either sim_X or sim_T are None the this function will attempt to plot the deterministic ODE
If we have 3 states or more, it will always be arrange such that it has 3 columns. Uses the operation from
odeutils
-
simulate_jump
(t, iteration, parallel=False, exact=False, full_output=False) Simulate the ode using stochastic simulation. It switches between a first reaction method and a \(\tau\)-leap algorithm internally. When a parallel backend exists, then a new random state (seed) will be used for each processor. This is due to a lack of appropriate parallel seed random number generator in python.
Parameters: - t: array like
the range of time points which we want to see the result of or the final time point
- iteration: int
number of iterations you wish to simulate
- parallel: bool, optional
Defaults to True
- exact: bool, optional
True if exact simulation is desired, defaults to False
- full_output: bool, optional
if we want additional information, sim_T
Returns: - sim_X: list
of length iteration each with (len(t),len(state)) if t is a vector, else it outputs unequal shape that was record of all the jumps
- sim_T: list or
numpy.ndarray
if t is a single value, it outputs unequal shape that was record of all the jumps. if t is a vector, it outputs t so that it is a
numpy.ndarray
instead
-
simulate_param
(t, iteration, parallel=False, full_output=False) Simulate the ode by generating new realization of the stochastic parameters and integrate the system deterministically.
Parameters: - t: array like
the range of time points which we want to see the result of
- iteration: int
number of iterations you wish to simulate
- parallel: bool, optional
Defaults to True
- full_output: bool, optional
if we want additional information, Y_all in the return, defaults to false
Returns: - Y:
numpy.ndarray
of shape (len(t), len(state)), mean of all the simulation
- Y_all:
np.ndarray
of shape (iteration, len(t), len(state))
-
total_transition
(state, t) Evaluate the total transition rate given state and time
Parameters: - state: array like
The current numerical value for the states which can be
numpy.ndarray
orlist
- t: double
The current time
Returns: - float
total rate
-
transition_matrix
(state, t) Evaluate the transition matrix given state and time
Parameters: - state: array like
The current numerical value for the states which can be
numpy.ndarray
orlist
- t: double
The current time
Returns: numpy.ndarray
a 2d array of size (M,M) where M is the number of transitions
-
transition_mean
(state, t) Evaluate the mean of the transitions given state and time. For m transitions and n states, we have
\[\begin{split}f_{j,k} &= \sum_{i=1}^{n} \frac{\partial a_{j}(x)}{\partial x_{i}} v_{i,k} \\ \mu_{j} &= \sum_{k=1}^{m} f_{j,k}(x)a_{k}(x) \\ \sigma^{2}_{j}(x) &= \sum_{k=1}^{m} f_{j,k}^{2}(x) a_{k}(x)\end{split}\]where \(v_{i,k}\) is the state change matrix.
Parameters: - state: array like
The current numerical value for the states which can be
numpy.ndarray
orlist
- t: double
The current time
Returns: numpy.ndarray
an array of size m where m is the number of transition
-
transition_var
(state, t) Evaluate the variance of the transitions given state and time
Parameters: - state: array like
The current numerical value for the states which can be
numpy.ndarray
orlist
- t: double
The current time
Returns: numpy.ndarray
an array of size M where M is the number of transition
-
transition_vector
(state, t) Evaluate the transition vector given state and time
Parameters: - state: array like
The current numerical value for the states which can be
numpy.ndarray
orlist
- t: double
The current time
Returns: numpy.ndarray
a 1d array of size K where K is the number of between states transitions and the number of birth death processes