toggleSwitch: The toggle switch model
In CharlotteJana/pdmpsim: Simulate Piecewise Deterministic Markov Processes

Description Usage Format Slots See Also Examples

This is an example for a piecewise deterministic markov process defined as pdmpModel. It models a gene regulation mechanism that is called toggle switch. This mechanism describes two genes A and B whose gene products are repressors to the other gene, i. e. the gene product of A blocks the gene expression of gene B and vice versa. This model is an example for a PDMP with two discrete and two continous variables. It is included to the package for demonstration purposes and is used in some unit tests and function examples.

1	toggleSwitch

An object of class pdmpModel.

init: There are two continous variables fA and fB which describe the concentration of gene products from gene A and gene B, respectivly. Both have initial value 0.5. The two discrete variables dA and dB describe the expression state of gene A and gene B and have initial value 1, which means that both genes are not blocked.
discStates: The discrete variables dA and dB both have codomain {0, 1} where 0 stands for "Gene A/B is blocked" and 1 stands for "Gene A/B is unblocked".
dynfunc: The dynamic of the continous variable fA depends on the state of discrete variable dA. In case dA = 0 it is given by dfA/dt = -bA⋅fA describing an exponential decay of the concentration of gene product A. If dA = 1, there is an additional term of linear growth leading to the ODE dfA/dt = -bA⋅fA + aA. Both formulas can be combined to dfA/dt = -bA⋅fA + aA⋅dA. Accordingly, the dynamic of dB is given as dfB/dt = -bB⋅fB + aB⋅dB.
jumpfunc: There are two jumptypes. If the jump is of type 1, only the value of dB is changed (from 0 to 1 or 1 to 0, respectivly). All other variables remain unchanged. The second jumptype changes only the value of dA.
ratefunc: A vector of length two determining the probability of a jump being of type 1 or 2. In this model, a molecule of gene product A can act as repressor and impede the gene expression of gene B. The rates for the first jumptype (which describes a change in the gene expression of gene B) are therefore k10A⋅fA for the change from the unblocked (dB = 1) to the blocked (dB = 0) state and k01A for the change from the blocked to the unblocked state of gene B. The rates for the second jumptype are generated in an analogous way because gene A is repressed by a molecule of gene product B and the second jumptype describes a change in the gene expression of gene A.
parms: There are are a number of parameters that appear in the dynamics and rates of the process, namely bA, bB, aA, aB, k01A, k10A, k01B and k10B. The values of the parameters in this example are artificial and are not based on real data.
times: The simulations will start at time t = 0 and end at t = 10 with step length 0.01.

simplePdmp for an easier example of a pdmpModel and pdmpModel for the formal description of the S4 class.

## the code used to generate this model:

toggleSwitch <- pdmpModel(
    descr = "Toggle Switch with two Promotors",
    parms = list(bA = 0.5, bB = 0.5, aA = 2, aB = 4, 
                 k01A = 0.5, k10A = 2, k01B = 0.3, k10B = 3),
    init = c(fA = 0.5, fB = 0.5, dA = 1.0, dB = 1.0),
    discStates = list(dA = c(0, 1), dB = c(0, 1)),
    times = c(from = 0, to = 100, by = 0.01),
    dynfunc = function(t, x, parms) {
       df <- with(as.list(c(x, parms)), c(-bA*fA + aA*dA, -bB*fB + aB*dB))
       return(c(df, 0, 0))
    },
    ratefunc = function(t, x, parms) {
       return(with(as.list(c(x, parms)), c(switch(dB+1, k01B, k10B*fA),
                                           switch(dA+1, k01A, k10A*fB))))
    },
    jumpfunc = function(t, x, parms, jtype){
       return(with(as.list(c(x, parms)), c(fA, fB, switch(jtype,
                                                          c(dA, 1-dB),
                                                          c(1-dA, dB)))))
    })

## load it and plot a simulation:
data("toggleSwitch")
plot(sim(toggleSwitch))