Skip to contents

Formula-based simulation constraints

Source: vignettes/run_until.Rmd
run_until.Rmd

Note: This article presents advanced topics on tissue simulation. Refer to vignette("tissue_simulation") for an introduction of the subject.

Disclaimer: RACES/rRACES internally implements the probability distributions using the C++11 random number distribution classes. The standard does not specify their algorithms, and the class implementations are left free for the compiler. Thus, the simulation output depends on the compiler used to compile RACES, and because of that, the results reported in this article may differ from those obtained by the reader.

rRACES implements a first order unquantified logic having variables representing the cardinality of the species, the number of events fired in a species (being duplications, deaths, or switches), and the simulation time. These variables and reals values are summed by +, subtracted by -, and multiplied by * and form expressions. The expressions are then compared with the standard semantics by >, >=, ==, !=, <=, and < to form relations. A formula in this language is either a relation, the conjunction of two formulas (&), or the disjunction of two formulas (|).

Any formula in above language expresses a condition on the simulation status and it can be used as the parameter of the method SpatialSimulation$run_until() to let the simulation evolve until the condition does not hold.

Variables

The variables represent one among the following quantities:

  • the cardinality of a species;
  • the number fired event among deaths, duplications and switches in a species;
  • the elapse simulation time.

All above variables can be built by using the method SpatialSimulation$var(). When the parameter is the string "Time", the elapsed simulation time variable is returned.

library(rRACES)

# set the seed of the random number generator
set.seed(0)

# build a spatial simulation and add two species to it
sim <- SpatialSimulation()
sim$add_mutant(name = "A",
               epigenetic_rates = c("+-" = 0.01, "-+" = 0.01),
               growth_rates = c("+" = 0.2, "-" = 0.08),
               death_rates = c("+" = 0.1, "-" = 0.01))

# get the variable representing the simulation time
v_time <- sim$var("Time")

v_time
#> Time

When the parameter is the name of a species, a variable representing the cardinality of the species is built.

# get the variable representing the cardinality of A+ in sim
va_p <- sim$var("A+")
va_p
#> |A+|

# get the variable representing the cardinality of A- in sim
va_m <- sim$var("A-")
va_m
#> |A-|

Finally, when the parameter is the name of a species followed by a . and the name of an event among deaths, duplications, or switches, SpatialSimulation$var() returns the variable associated with the number of the corresponding event in the specified species.

# get the variable representing the number of epigenetic
# switches from A+
va_ps <- sim$var("A+.switches")
va_ps
#> |A+.switches|

# get the variable representing the number of duplications in A+
sim$var("A+.duplications")
#> |A+.duplications|

# get the variable representing the number of deaths in A+
sim$var("A+.deaths")
#> |A+.deaths|

Expressions and Formulas

An expression is one of the following object:

  • a variable, e.g., sim$var("A+");
  • a numeric value, e.g., 3.4;
  • the sum of two expressions, e.g., sim$var("A+") + 3.4;
  • the subtraction of two expressions, e.g., sim$var("A+") - 3.4;
  • the multiplication of two expressions, e.g., sim$var("A+") * 3.4.

Two expression can be related by <=, <, ==, !=, > and >=.

A formula is:

  • a relation among two expressions, e.g., sim$var("A+")>=2;
  • the conjunction of two formulas, e.g., sim$var("A+")>=2 & sim$var("A+")<=500;
  • the disjunction of two formulas, e.g., sim$var("A+")>=2 | sim$var("A+")<=500.

The method SpatialSimulation$run_until()

The method SpatialSimulation$run_until() takes as the parameter a formula and lets the simulation evolve until the formula does not hold.

# build a condition stating that the cardinality of A+ doubles
# that of A-
c1 <- va_p >= 2 * va_m
c1
#> |A+|>=2*|A-|

# build a condition that holds when there are more than
# 100000 live cells of mutant A
c2 <- va_p + va_m > 1e5
c2
#> |A+|+|A-|>100000

# build a condition that holds when less than 4000
# epigenetic switches from the species A+ have occured
c3 <- va_ps < 4000
c3
#> 4000>|A+.switches|

# build a condition that holds when 40 time unit have been
# simulated at least
c4 <- v_time >= 40
c4
#> Time>=40

# build a condition that holds when c4 and at least one
# among c1, c2, and c3 hold
c5 <- c4 & (c1 | c2 | c3)
c5
#> Time>=40 and (|A+|>=2*|A-| or |A+|+|A-|>100000 or 4000>|A+.switches|)

# place the initial cell
sim$place_cell("A+", 500, 500)

# run the simulation while c5 does not hold
sim$run_until(c5)
#>  [████████████████████████████████████████] 100% [00m:00s] Saving snapshot

sim
#> ──  rRACES   D   S   M  /var/folders/tb/jqmdpgxs2t5129bny6pb96680000gn/T/rRACES_
#> 
#> ── Species: 2, with epigenetics
#>    
#>    =======  ====  ====  ====  ======  ===
#>    species   λ      δ    ε    counts   % 
#>    =======  ====  ====  ====  ======  ===
#>         A-  0.08  0.01  0.01    30    75 
#>         A+  0.20  0.10  0.01    10    25 
#>    =======  ====  ====  ====  ======  ===
#> 
#> ── Firings: 73 total
#> 
#>  Species [A-]:  5 (deaths), 37 (duplications) and  5 (switches)
#>  Species [A+]:  8 (deaths), 15 (duplications) and  3 (switches)
#> ✖ The simulation has no samples yet!
sim$get_clock()
#> [1] 40.34598