This R package provides an interface to define compartmental infectious disease models, typeset the model in LaTeX, simulate the model numerically by solving systems of ODEs, or simulate the model stochastically using the Gillespie method. If you are simply interested in simulating an ODE model, this package also provides a convenient interface to define an ODE model with ease.
REpiSim supports two types of models:
- Generic ODE models, defined directly by differential equations.
- Compartmental models, defined by transitions between states.
The R6 class used to define a generic ODE model is Model.
Equations may be defined either:
- directly in the
new()method, or - added one at a time using the
compartment()method.
Parameters appearing in the equations are automatically detected. Algebraic substitutions can be defined using the where() method.
An equation is specified using R formula syntax. For example,
S ~ -beta*S*Idefines a compartment S with derivative
[ \frac{dS}{dt} = -\beta S I. ]
A basic SIR model can be defined as follows:
m <- Model$new(
S ~ -beta*S*I/N,
I ~ beta*S*I/N - gamma*I,
R ~ gamma*I
)
print(m)Here N appears as a parameter. It can instead be defined as a substitution:
m$where(N = S + I + R)
print(m)Differential equations and substitutions can be typeset using a TexFormatter object:
tex <- TexFormatter$new()
tex$typeset.equation(m$equations$equations)
tex$typeset.equation(m$equations$where)Individual expressions can also be typeset:
tex$typeset(quote(beta*S*I))Custom LaTeX symbols can be defined using set.symbols().
Numerical simulation of an ODE model is handled by the ODE class, which wraps the deSolve package.
sim <- ODE$new(m)
x <- sim$simulate(
t = 0:100,
y0 = c(S = 1000, I = 10, R = 0),
parms = c(beta = 0.4, gamma = 0.2)
)
head(x)The returned object is a data frame containing the simulated trajectories. You may select a subset of variables via the vars argument, and substitution values can also be returned.
For many infectious disease models, it is more natural to describe the dynamics in terms of state transitions (for example, Susceptible → Infected → Recovered) rather than writing differential equations explicitly.
REpiSim provides the Compartmental class for this purpose. The Compartmental class is a subclass of Model. Internally, all transition rules are automatically translated into a system of ODEs, so everything that works for Model (LaTeX typesetting, numerical simulation, etc.) also works for Compartmental.
Compartments (state variables) are declared when creating the model:
m <- Compartmental$new(S, I, R)This declares S, I, and R as state variables, with their derivatives initially set to zero.
Transitions are added using the transition() method. A transition has the general form:
from -> to ~ rateA standard SIR model can be written as:
m$transition(S -> I ~ beta*S*I/N, N = S + I + R)
m$transition(I -> R ~ gamma*I)This corresponds to the ODE system
[ \begin{aligned} \frac{dS}{dt} &= -\beta S I / N, \ \frac{dI}{dt} &= \beta S I / N - \gamma I, \ \frac{dR}{dt} &= \gamma I. \end{aligned} ]
Model parameters such as beta and gamma are automatically detected. Substitutions (for example N = S + I + R) may be supplied directly in the transition call.
Births, deaths, or other external flows can be modeled using NULL as the source or destination:
- Births into a compartment:
m$transition(NULL -> S ~ mu)
- Deaths from a compartment:
m$transition(I -> NULL ~ mu*I)
For per-capita transition rates, you may either write them explicitly or use the percapita option:
m$transition(I -> R ~ gamma, percapita = TRUE)This is interpreted as a transition rate equal to gamma * I.
Because Compartmental is a subclass of Model, the generated differential equations can be inspected directly:
print(m)They can also be typeset using TexFormatter.
Once defined, compartmental models can be simulated in the same way as generic models.
Deterministic (ODE) simulation:
sim <- ODE$new(m)
x <- sim$simulate(
t = 0:100,
y0 = c(S = 1000, I = 10, R = 0),
parms = c(beta = 0.4, gamma = 0.2)
)Stochastic (Gillespie) simulation:
sim <- RGillespie$new(m)
x <- sim$simulate(
t = 0:100,
y0 = c(S = 1000, I = 10, R = 0),
parms = c(beta = 0.4, gamma = 0.2)
)- Use
Modelwhen you want full control over the differential equations or are working with non-compartmental systems. - Use
Compartmentalwhen your model is naturally described by transitions between states; this approach is often more concise and less error-prone.
REpiSim can generate TikZ flowcharts for compartmental models using the flowchart() function.
REpiSim provides tools to:
- define ODE models using symbolic equations,
- define compartmental models using transition rules,
- typeset equations in LaTeX,
- simulate models deterministically and stochastically, and
- generate flowcharts for visualization.