Systems

Example Systems included in PathWeightSampling.jl

PathWeightSampling.gene_expression_systemFunction
gene_expression_system(;
    mean_s=50,
    mean_x=mean_s,
    corr_time_s=1.0,
    corr_time_x=0.1,
    kappa=mean_s / corr_time_s,
    lambda=1 / corr_time_s,
    mu=1 / corr_time_x,
    rho=mu * mean_x / mean_s,
    u0=[mean_s, mean_x],
    dtimes=0:0.1:2.0
)

Creates a system for a very simple model of gene expression.

Model Description

This model consists of four reactions:

                κ
reaction 1:  ∅ ---> S

                λ
reaction 2:  S ---> ∅

                ρ
reaction 3:  S ---> S + X

                μ
reaction 4:  X ---> ∅

The first two reactions specify the evolution of the input signal S, and last two reactions specify the evolution of the output X. Thus, both the input and output signal are modeled as a simple birth-death process, however the birth rate of X increases with higher copy numbers of S.

Examples

The values of the reaction rates can be specified directly as follows:

import PathWeightSampling as PWS
system = PWS.gene_expression_system(kappa = 10.0, lambda = 0.1, rho = 1.0, mu = 1.0)

# output

MarkovJumpSystem with 2 species and 4 reactions
ReactionSet with species S, X
k1 =  10.00: ∅ ---> S
k2 =   0.10: S ---> ∅
k3 =   1.00: S ---> S + X
k4 =   1.00: X ---> ∅

Initial condition:
    S = 50
    X = 50

Alternatively, the reaction rates can be specified indirectly through the following arguments:

  • mean_s: The average copy number of S when the system has relaxed to steady state.
  • mean_x: The average copy number of X when the system has relaxed to steady state.
  • corr_time_s: The input signal correlation time. The shorter this time is, the faster the fluctuations in the input signal.
  • corr_time_x: The output signal correlation time. This sets the timescale of output fluctuations.
import PathWeightSampling as PWS
system = PWS.gene_expression_system(mean_s=25, mean_x=50, corr_time_s=1.0, corr_time_x=0.75)

# output

MarkovJumpSystem with 2 species and 4 reactions
ReactionSet with species S, X
k1 =  25.00: ∅ ---> S
k2 =   1.00: S ---> ∅
k3 =   2.67: S ---> S + X
k4 =   1.33: X ---> ∅

Initial condition:
    S = 25
    X = 50
source
PathWeightSampling.chemotaxis_systemFunction
chemotaxis_system(;
    n_clusters=25,
    n_chey=10000,
    methylation_sites=4,
    duration=200.0,
    dt=0.1,
    sde_dt=dt / 10,
    c_0=100.0,
    Kₐ=2900.0,
    Kᵢ=18.0,
    n=15, # cooperativity
    k_R=0.1,
    k_B=0.2,
    a_0=k_R / (k_R + k_B),
    δf=-2.0,
    m_0=0.5 * n,
    k_Z=10.0,
    phi_y=1 / 6,
    k_A=k_Z * phi_y / ((1 - phi_y) * a_0 * n_clusters),
    velocity_decay=0.862,
    velocity_noise=sqrt(2 * velocity_decay * 157.1),
    gradient_steepness=0.2e-3,
)

Create a system for a stochastic chemotaxis model.

Model description

This model describes the bacterial chemotaxis signaling network.

Rml -> Rm(l+1), with rate kon L(t) (lmax - l): ligand binding to active state
Rm(l+1) -> Rml, with rate koff_A l: ligand unbinding
Rml -> R(m+1)l with rate k_R: methylation rate
R(m+1)l -> Rml with rate k_B: demethylation rate

Example

import PathWeightSampling as PWS
system = PWS.chemotaxis_system()

# output

HybridJumpSystem with 64 species and 184 reactions
PathWeightSampling.ChemotaxisJumps(2900.0, 18.0, 15, 0.2, 0.1, 0.23999999999999996, 10.0, -2.0, 7.5, 1, 2:62, 63, 64)

Initial condition:
    L = 100.0
    R0 = 0.0
    R1 = 0.0
    R2 = 0.0
    R3 = 0.0
    R4 = 0.0
    R5 = 0.0
    R6 = 0.0
    R7 = 0.0
    R8 = 0.0
    R9 = 0.0
    R10 = 0.0
    R11 = 0.0
    R12 = 0.0
    R13 = 0.0
    R14 = 0.0
    R15 = 0.0
    R16 = 0.0
    R17 = 0.0
    R18 = 0.0
    R19 = 1.0
    R20 = 8.0
    R21 = 11.0
    R22 = 5.0
    R23 = 0.0
    R24 = 0.0
    R25 = 0.0
    R26 = 0.0
    R27 = 0.0
    R28 = 0.0
    R29 = 0.0
    R30 = 0.0
    R31 = 0.0
    R32 = 0.0
    R33 = 0.0
    R34 = 0.0
    R35 = 0.0
    R36 = 0.0
    R37 = 0.0
    R38 = 0.0
    R39 = 0.0
    R40 = 0.0
    R41 = 0.0
    R42 = 0.0
    R43 = 0.0
    R44 = 0.0
    R45 = 0.0
    R46 = 0.0
    R47 = 0.0
    R48 = 0.0
    R49 = 0.0
    R50 = 0.0
    R51 = 0.0
    R52 = 0.0
    R53 = 0.0
    R54 = 0.0
    R55 = 0.0
    R56 = 0.0
    R57 = 0.0
    R58 = 0.0
    R59 = 0.0
    R60 = 0.0
    Yp = 0.0
    Y = 10000.0
source