Package 'SimInf'

Description

Approximate Bayesian computation

Usage

abc(
  model,
  priors = NULL,
  n_particles = NULL,
  n_init = NULL,
  distance = NULL,
  tolerance = NULL,
  data = NULL,
  verbose = FALSE,
  post_gen = NULL,
  init_model = NULL
)

## S4 method for signature 'SimInf_model'
abc(
  model,
  priors = NULL,
  n_particles = NULL,
  n_init = NULL,
  distance = NULL,
  tolerance = NULL,
  data = NULL,
  verbose = FALSE,
  post_gen = NULL,
  init_model = NULL
)

Arguments

Value

A SimInf_abc object.

References

T. Toni, D. Welch, N. Strelkowa, A. Ipsen, and M. P. H. Stumpf. Approximate Bayesian computation scheme for parameter inference and model selection in dynamical systems. Journal of the Royal Society Interface 6, 187–202, 2009. doi:10.1098/rsif.2008.0172

U. Simola, J. Cisewski-Kehe, M. U. Gutmann, J. Corander. Adaptive Approximate Bayesian Computation Tolerance Selection. Bayesian Analysis, 16(2), 397–423, 2021. doi:10.1214/20-BA1211

See Also

continue_abc.

Examples

## Not run: 
## Let us consider an SIR model in a closed population with N = 100
## individuals of whom one is initially infectious and the rest are
## susceptible. First, generate one realisation (with a specified
## seed) from the model with known parameters \code{beta = 0.16} and
## \code{gamma = 0.077}. Then, use \code{abc} to infer the (known)
## parameters from the simulated data.
model <- SIR(u0 = data.frame(S = 99, I = 1, R = 0),
             tspan = 1:100,
             beta = 0.16,
             gamma = 0.077)

## Run the SIR model and plot the number of infectious.
set.seed(22)
infectious <- trajectory(run(model), "I")$I
plot(infectious, type = "s")

## The distance function to accept or reject a proposal. Each node
## in the simulated trajectory (contained in the 'result' object)
## represents one proposal.
distance <- function(result, ...) {
    ## Extract the time-series of infectious in each node as a
    ## data.frame.
    sim <- trajectory(result, "I")

    ## Split the 'sim' data.frame by node and calculate the sum of the
    ## squared distance at each time-point for each node.
    dist <- tapply(sim$I, sim$node, function(sim_infectious) {
        sum((infectious - sim_infectious)^2)
    })

    ## Return the distance for each node. Each proposal will be
    ## accepted or rejected depending on if the distance is less than
    ## the tolerance for the current generation.
    dist
}

## Fit the model parameters using ABC-SMC and adaptive tolerance
## selection. The priors for the parameters are specified using a
## formula notation. Here we use a uniform distribtion for each
## parameter with lower bound = 0 and upper bound = 1. Note that we
## use a low number particles here to keep the run-time of the example
## short. In practice you would want to use many more to ensure better
## approximations.
fit <- abc(model = model,
           priors = c(beta ~ uniform(0, 1), gamma ~ uniform(0, 1)),
           n_particles = 100,
           n_init = 1000,
           distance = distance,
           verbose = TRUE)

## Print a brief summary.
fit

## Display the ABC posterior distribution.
plot(fit)

## End(Not run)

Description

A utility function to augment local model parameters (ldata) with spatial coupling data from neighboring nodes.

Usage

add_spatial_coupling_to_ldata(
  x,
  y,
  cutoff,
  ldata = NULL,
  min_dist = NULL,
  na_fail = TRUE
)

Arguments

Details

The function calculates distances between nodes based on projected coordinates (x, y) and appends neighbor data to the ldata matrix for each node.

Output Format: The returned matrix has the same number of rows as the input ldata. The columns are organized as follows:

• Local Parameters: The first $n$ columns correspond to the original local parameters passed in ldata.

• Neighbor Pairs: Following the local parameters, the data is stored as pairs of columns: (neighbor_index, distance). The neighbor_index is a zero-based index of the neighbor node. The distance is the Euclidean distance to that neighbor.

• Stop Marker: Each node's neighbor list is terminated by a pair (-1, 0) in the position where the next neighbor_index would appear. This marker appears exactly once per node, immediately after the last neighbor.

Value

A numeric matrix with the same number of rows as ldata, but with additional columns containing the neighbor indices and distances as described in the ‘Output Format’ section.

See Also

distance_matrix for computing pairwise distances between nodes, and edge_properties_to_matrix for a similar utility that converts edge properties to a matrix format.

Examples

## Generate a 5 x 5 grid of nodes separated by 1000m.
nodes <- expand.grid(
  x = seq(from = 0, to = 4000, by = 1000),
  y = seq(from = 0, to = 4000, by = 1000)
)

## Create local data with one parameter per node.
ldata <- matrix(0.1, nrow = 1, ncol = nrow(nodes))

## Add spatial coupling with a 2500m cutoff.
ldata_augmented <- add_spatial_coupling_to_ldata(
  x = nodes$x,
  y = nodes$y,
  cutoff = 2500,
  ldata = ldata
)

## Inspect the result for the first node.
ldata_augmented[, 1]

Description

Convert the results of an Approximate Bayesian Computation (ABC) analysis into a single data.frame. This function extracts the particle parameters, their acceptance weights, and the generation number for every particle across all generations.

Usage

## S3 method for class 'SimInf_abc'
as.data.frame(x, ...)

Arguments

Details

The resulting data.frame has one row per particle. The columns include:

• generation: The generation number (integer).

• weight: The normalized weight of the particle (numeric).

• Parameter columns: One column for each parameter estimated in the ABC analysis (e.g., beta, gamma, sigma). The column names match the parameter names defined in the priors.

Value

A data.frame containing all particles from all generations.

See Also

abc for running the ABC analysis, SimInf_abc for the class definition, and continue_abc for continuing an existing ABC run.

Description

Convert the scheduled events stored in a SimInf_events object into a data.frame. This function extracts the event type, time, source node, destination node, number of individuals, proportion, and the specific columns from the select (E) and shift (N) matrices that define how each event modifies the compartment state. The resulting data.frame has one row per scheduled event.

Usage

## S3 method for class 'SimInf_events'
as.data.frame(x, ...)

Arguments

Value

A data.frame with columns:

• event: Event type (integer or character, depending on input).

• time: Time of the event (integer or Date, depending on input).

• node: Source node identifier.

• dest: Destination node identifier (may be NA).

• n: Number of individuals affected.

• proportion: Proportion of the population affected (if applicable).

• select: The column vector from the select matrix (E) that defines how the event modifies the compartment state.

• shift: The column vector from the shift matrix (N) that defines how the event modifies the compartment state.

See Also

SimInf_events for the class definition and events for extracting events from a model.

Examples

## Create an 'SIR' model with 1600 cattle herds (nodes) and
## initialize it to run over 4*365 days. Define 'tspan' to record
## the state of the system at daily time-points. Load scheduled
## events for the population of nodes with births, deaths and
## between-node movements of individuals.
model <- SIR(
  u0     = u0_SIR(),
  tspan  = seq(from = 1, to = 4*365, by = 1),
  events = events_SIR(),
  beta   = 0.16,
  gamma  = 0.01
)

## Extract the events from the model and convert to a data frame.
head(as.data.frame(events(model)))

Description

Convert the cleaned individual-level events stored in a SimInf_individual_events object into a data.frame. This function extracts the individual identifier, event type, time, source node, and destination node. The resulting data.frame has one row per event.

Usage

## S3 method for class 'SimInf_individual_events'
as.data.frame(x, ...)

Arguments

Value

A data.frame with columns:

• id: Identifier of the individual (integer or character, depending on input).

• event: Event type (integer or character, depending on input).

• time: Time of the event (numeric or Date, depending on input).

• node: Source node identifier (integer or character, depending on input).

• dest: Destination node identifier (integer or character, depending on input) (may be NA).

See Also

SimInf_individual_events for the class definition and individual_events for cleaning livestock event data and prepare it for usage in SimInf.

Description

Extract the posterior samples from the MCMC chain stored in a SimInf_pmcmc object and convert them into a data.frame.

Usage

## S3 method for class 'SimInf_pmcmc'
as.data.frame(x, ...)

Arguments

Details

The resulting data.frame contains one row per MCMC iteration and one column per parameter. These samples represent the joint posterior distribution of the parameters. This format is convenient for post-processing and visualization.

Value

A data.frame where rows represent MCMC iterations and columns represent the posterior samples of the parameters.

See Also

pmcmc for running the PMCMC analysis, SimInf_pmcmc for the class definition, and continue_pmcmc for continuing an existing PMCMC run.

Description

Produce box-and-whisker plots summarizing the distribution of the number of individuals in specified model compartments. The plots aggregate data across all time points in tspan and the selected nodes (specified by index).

Usage

## S4 method for signature 'SimInf_model'
boxplot(x, compartments = NULL, index = NULL, ...)

Arguments

Details

This function is useful for visualizing the variability and range of compartment counts across the entire simulation trajectory.

Examples

## For reproducibility, set the seed and number of threads.
set.seed(123)
set_num_threads(1)

## Create an 'SIR' model with 10 nodes.
model <- SIR(
  u0 = data.frame(
    S = rep(99, 10),
    I = rep(1, 10),
    R = rep(0, 10)
  ),
  tspan = 1:100,
  beta = 0.16,
  gamma = 0.077
)

## Run the model.
result <- run(model)

## Create a boxplot for all compartments across all nodes and time
## points.
boxplot(result)

## Create a boxplot for the S and I compartments, but only for
## nodes 1 and 2.
boxplot(result, compartments = c("S", "I"), index = 1:2)

Description

This function extracts the C source code associated with the model. The code may have been automatically generated by mparse or manually provided by the user during model creation. It is useful for inspecting the implementation, debugging complex propensity functions, or verifying the code used in the simulation.

Usage

C_code(model)

Arguments

Value

A character vector where each element corresponds to a line of the C source code stored in the model.

Examples

## Extract code from an mparse-generated SIR model.
## Using writeLines formats the character vector output for
## readability.
model <- mparse(
  transitions = c("S -> beta * S * I/(S + I + R) -> I",
                  "I -> gamma * I -> R"),
  compartments = c("S", "I", "R"),
  gdata = c(beta = 0.16, gamma = 0.077),
  u0 = data.frame(S = 99, I = 1, R = 0),
  tspan = 1:10,
  use_enum = TRUE
)

writeLines(C_code(model))

Description

Run more generations of ABC SMC

Usage

continue_abc(
  object,
  tolerance = NULL,
  data = NULL,
  verbose = FALSE,
  post_gen = NULL
)

## S4 method for signature 'SimInf_abc'
continue_abc(
  object,
  tolerance = NULL,
  data = NULL,
  verbose = FALSE,
  post_gen = NULL
)

Arguments

Value

A SimInf_abc object.

Description

Extends a previously computed PMCMC chain by running additional iterations of the Metropolis-Hastings sampler, continuing from the final state of the previous chain.

Usage

continue_pmcmc(
  object,
  obs_process,
  n_iterations,
  post_proposal = NULL,
  init_model = NULL,
  post_particle = NULL,
  verbose = FALSE
)

## S4 method for signature 'SimInf_pmcmc'
continue_pmcmc(
  object,
  obs_process,
  n_iterations,
  post_proposal = NULL,
  init_model = NULL,
  post_particle = NULL,
  verbose = FALSE
)

Arguments

Value

The updated SimInf_pmcmc object with the chain extended by n_iterations new rows.

See Also

pmcmc for initiating a new PMCMC chain.

Description

Calculate the Euclidean distances between all pairs of nodes based on their projected coordinates (x, y). Distances greater than the specified cutoff are excluded from the result (stored as zeros in the sparse matrix).

Usage

distance_matrix(x, y, cutoff, min_dist = NULL, na_fail = TRUE)

Arguments

Details

The result is a symmetric sparse matrix (dgCMatrix) where the element d[i, j] contains the distance between node i and node j if it is less than or equal to cutoff, and 0 otherwise.

Value

A symmetric sparse matrix of class dgCMatrix. Non-zero entries represent distances $\le$ cutoff; entries outside the cutoff are implicitly zero.

Examples

## Generate a 10 x 10 grid of nodes separated by 100m.
nodes <- expand.grid(
  x = seq(from = 0, to = 900, by = 100),
  y = seq(from = 0, to = 900, by = 100)
)
plot(nodes, main = "Node Grid", asp = 1)

## Calculate distances with a 300m cutoff.
## Only neighbors within 300m will have non-zero entries.
d <- distance_matrix(
  x = nodes$x,
  y = nodes$y,
  cutoff = 300
)

## Inspect the sparse matrix structure.
## Note: The matrix is symmetric and the diagonal is zero.
d[1:10, 1:10]

## Count the number of neighbors for the first node.
sum(d[1, ] > 0)

Description

A utility function to convert a data frame of edge properties (e.g., movement rates, counts) into the specific matrix format required for ldata in spatial models. This format allows the C code to efficiently iterate over all incoming edges for a specific node.

Usage

edge_properties_to_matrix(edges, n_nodes)

Arguments

Details

The function converts the edges data frame into a numeric matrix where:

• Each column corresponds to a to node.

• Each column contains a single sequence of data blocks, one for each incoming edge to that node.

• Each block starts with the zero-based index of the from node.

• The subsequent rows in the block contain the property values (e.g., rate, count).

• The sequence for a node ends with a stop marker (-1) in the first row of the block (the position where the from index is stored). This marker appears exactly once per column, immediately after the last incoming edge.

• Unused cells in the matrix (below the stop marker) are filled with NaN.

The edges data frame must contain columns from and to with valid node indices (1-based). All edges pointing to the same to node must be unique (no duplicate from -> to pairs).

Value

A numeric matrix with dimensions determined by the maximum number of incoming edges for any node. Columns correspond to to nodes. Entries are NaN where no data exists.

Examples

## Define edge properties: from, to, rate, and count.
edges <- data.frame(
  from  = c(2, 3, 4, 1, 4, 5, 1, 3, 1, 3),
  to    = c(1, 1, 1, 2, 3, 3, 4, 4, 5, 5),
  rate  = c(0.2, 0.01, 0.79, 1, 0.2, 0.05, 0.2, 0.8, 0.2, 0.8),
  count = c(5, 5, 5, 50, 10, 10, 5, 5, 5, 5)
)

## Convert to matrix for 6 nodes.
mat <- edge_properties_to_matrix(edges, 6)
print(mat)

## Interpretation of Column 1 (to node 1): The column contains a
## single sequence of 3 blocks (edges from # 2, 3, and 4).
## Row 1: Index of 'from' node 2 (2-1 = 1).
## Row 2: rate=0.2
## Row 3: count=5
## Row 4: Index of 'from' node 3 (3-1 = 2).
## Row 5: rate=0.01
## Row 6: count=5
## Row 7: Index of 'from' node 4 (4-1 = 3).
## Row 8: rate=0.79
## Row 9: count=5
## Row 10: Stop marker (-1) indicating the end of the list for node 1.

Description

Retrieve the SimInf_events object containing the schedule of discrete events (e.g., births, deaths, movements) associated with a SimInf_model. This object holds the timing, location, and type of each event, as well as the matrices defining how events affect the model state.

Usage

events(object, ...)

## S4 method for signature 'SimInf_model'
events(object, ...)

Arguments

Value

A SimInf_events object containing the event schedule and associated matrices (E and N).

See Also

SimInf_events for details on the structure of the returned event object (slots E (select matrix), N (shift matrix), event, etc.). events_SIR, events_SEIR, events_SISe3 for examples of pre-defined event datasets. mparse for defining custom models with event schedules. run for executing the simulation with the scheduled events. Vignette "Scheduled events" for a comprehensive tutorial on defining event data, using the select (E) and shift (N) matrices, and simulating complex demographic and movement processes.

Examples

## Create an SIR model with scheduled events.
model <- SIR(
  u0     = u0_SIR(),
  tspan  = 1:(4 * 365),
  events = events_SIR(),
  beta   = 0.16,
  gamma  = 0.077
)

## Extract the events and display a summary.
ev <- events(model)
summary(ev)

## Plot the event schedule over time.
plot(ev)

Description

Dataset containing 466,692 scheduled events for a population of 1,600 cattle herds over 1,460 days (4 years). Demonstrates how demographic and movement events affect SEIR dynamics in a cattle disease context.

Usage

events_SEIR()

Details

The event data contains three types of scheduled events that affect cattle herds (nodes):

Exit

Deaths or removal of cattle from a herd (n = 182,535). These events decrease the population and affect all disease compartments proportionally.

Enter

Births or introduction of cattle to a herd (n = 182,685). These events add susceptible cattle to herds, increasing overall herd size.

External transfer

Movement of cattle between herds (n = 101,472). These events transfer cattle from one herd to another, potentially facilitating between-herd disease transmission.

The select column in the returned data frame is mapped to the columns of the internal select matrix:

• select = 1 corresponds to Enter events, targeting the Susceptible (S) compartment.

• select = 2 corresponds to Exit and External Transfer events, targeting all compartments.

Events are distributed across all 1,600 herds over the 4-year period. These are synthetic data generated to illustrate how to incorporate scheduled events (such as births, deaths, and movements) into a compartment model in the SimInf framework.

Value

A data.frame with columns:

event

Event type: "exit", "enter", or "extTrans".

time

Day when event occurs (1-1460).

node

Affected herd identifier (1-1600).

dest

Destination herd for external transfer events.

n

Number of cattle affected.

proportion

0. Not used in this example.

select

Model compartment to affect (see SimInf_events).

shift

0. Not used in this example.

See Also

u0_SEIR for the corresponding initial cattle population, SEIR for creating SEIR models with these events, and SimInf_events for event structure details

Examples

## For reproducibility, call the set.seed() function and specify the
## number of threads to use. To use all available threads, remove the
## set_num_threads() call.
set.seed(123)
set_num_threads(1)

## Create a 'SEIR' model with 1600 cattle herds (nodes) and initialize
## it to run over 4*365 days. Add ten exposed animals to the first
## herd. Define 'tspan' to record the state of the system at weekly
## time-points. Load scheduled events for the population of nodes with
## births, deaths and between-node movements of individuals.
u0 <- u0_SEIR()
u0$E[1] <- 10
model <- SEIR(
    u0      = u0,
    tspan   = seq(from = 1, to = 4*365, by = 7),
    events  = events_SEIR(),
    beta    = 0.16,
    epsilon = 0.25,
    gamma   = 0.01
)

## Display the number of cattle affected by each event type per day.
plot(events(model))

## Run the model to generate a single stochastic trajectory.
result <- run(model)

## Plot the median and interquartile range of the number of
## susceptible, exposed, infected and recovered individuals.
plot(result)

## Plot the trajectory for the first herd.
plot(result, index = 1)

## Summarize the trajectory. The summary includes the number of events
## by event type.
summary(result)

Description

Dataset containing 466,692 scheduled events for a population of 1,600 cattle herds over 1,460 days (4 years). Demonstrates how demographic and movement events affect SIR dynamics in a cattle disease context.

Usage

events_SIR()

Details

The event data contains three types of scheduled events that affect cattle herds (nodes):

Exit

Deaths or removal of cattle from a herd (n = 182,535). These events decrease the population and remove cattle from the disease system.

Enter

Births or introduction of cattle to a herd (n = 182,685). These events add susceptible cattle to herds, increasing potential targets for infection.

External transfer

Movement of cattle between herds (n = 101,472). These events transfer cattle from one herd to another, potentially spreading disease across the herd network.

The select column in the returned data frame is mapped to the columns of the internal select matrix (select_matrix_SIR):

• select = 1 corresponds to Enter events, targeting the Susceptible (S) compartment.

• select = 4 corresponds to Exit and External Transfer events, targeting all compartments (S, I, and R).

Events are distributed across all 1,600 herds over the 4-year period. These are synthetic data generated to illustrate how to incorporate scheduled events (such as births, deaths, and movements) into a compartment model in the SimInf framework.

Value

A data.frame with columns:

event

Event type: "exit", "enter", or "extTrans".

time

Day when event occurs (1-1460).

node

Affected herd identifier (1-1600).

dest

Destination herd for external transfer events.

n

Number of cattle affected.

proportion

0. Not used in this example.

select

Model compartment to affect (see SimInf_events).

shift

0. Not used in this example.

See Also

u0_SIR for the corresponding initial cattle population, SIR for creating SIR models with these events, and SimInf_events for event structure details

Examples

## For reproducibility, call the set.seed() function and specify the
## number of threads to use. To use all available threads, remove the
## set_num_threads() call.
set.seed(123)
set_num_threads(1)

## Create an 'SIR' model with 1600 cattle herds (nodes) and initialize
## it to run over 4*365 days. Add one infected animal to the first
## herd to seed the outbreak. Define 'tspan' to record the state of
## the system at daily time-points. Load scheduled events for the
## population of nodes with births, deaths and between-node movements
## of individuals.
u0 <- u0_SIR()
u0$I[1] <- 1
model <- SIR(
    u0     = u0,
    tspan  = seq(from = 1, to = 4*365, by = 1),
    events = events_SIR(),
    beta   = 0.16,
    gamma  = 0.01
)

## Display the number of cattle affected by each event type per day.
plot(events(model))

## Run the model to generate a single stochastic trajectory.
result <- run(model)

## Plot the median and interquartile range of the number of
## susceptible, infected and recovered individuals.
plot(result)

## Plot the trajectory for the first herd.
plot(result, index = 1)

## Summarize the trajectory. The summary includes the number of events
## by event type.
summary(result)

Description

Dataset containing 466,692 scheduled events for a population of 1,600 cattle herds over 1,460 days (4 years). Demonstrates how demographic and movement events affect SIS dynamics in a cattle disease context.

Usage

events_SIS()

Details

The event data contains three types of scheduled events that affect cattle herds (nodes):

Exit

Deaths or removal of cattle from a herd (n = 182,535). These events decrease the population in both susceptible and infected compartments.

Enter

Births or introduction of cattle to a herd (n = 182,685). These events add susceptible cattle to herds.

External transfer

Movement of cattle between herds (n = 101,472). These events transfer cattle from one herd to another, potentially spreading disease across the herd network. Either susceptible or infected animals may be transferred.

The select column in the returned data frame is mapped to the columns of the internal select matrix:

• select = 1 corresponds to Enter events, targeting the Susceptible (S) compartment.

• select = 2 corresponds to Exit and External Transfer events, targeting all compartments (S and I).

Events are distributed across all 1,600 herds over the 4-year period. These are synthetic data generated to illustrate how to incorporate scheduled events (such as births, deaths, and movements) into a compartment model in the SimInf framework.

Value

A data.frame with columns:

event

Event type: "exit", "enter", or "extTrans".

time

Day when event occurs (1-1460).

node

Affected herd identifier (1-1600).

dest

Destination herd for external transfer events.

n

Number of cattle affected.

proportion

0. Not used in this example.

select

Model compartment to affect (see SimInf_events).

shift

0. Not used in this example.

See Also

u0_SIS for the corresponding initial cattle population, SIS for creating SIS models with these events, and SimInf_events for event structure details

Examples

## For reproducibility, call the set.seed() function and specify the
## number of threads to use. To use all available threads, remove the
## set_num_threads() call.
set.seed(123)
set_num_threads(1)

## Create an 'SIS' model with 1600 cattle herds (nodes) and initialize
## it to run over 4*365 days. Add one infected animal to the first
## herd to seed the outbreak. Define 'tspan' to record the state of
## the system at daily time-points. Load scheduled events for the
## population of nodes with births, deaths and between-node movements
## of individuals.
u0 <- u0_SIS()
u0$I[1] <- 1
model <- SIS(
    u0     = u0,
    tspan  = seq(from = 1, to = 4*365, by = 1),
    events = events_SIS(),
    beta   = 0.16,
    gamma  = 0.01
)

## Display the number of cattle affected by each event type per day.
plot(events(model))

## Run the model to generate a single stochastic trajectory.
result <- run(model)

## Plot the median and interquartile range of the number of
## susceptible and infected individuals.
plot(result)

## Plot the trajectory for the first herd.
plot(result, index = 1)

## Summarize the trajectory. The summary includes the number of events
## by event type.
summary(result)

Description

Dataset containing 783,773 scheduled events for a population of 1,600 cattle herds stratified by age over 1,460 days (4 years). Demonstrates how demographic, movement, and age-transition events affect SISe3 dynamics in a cattle disease context.

Usage

data(events_SISe3)

Format

A data.frame

Details

This dataset contains four types of scheduled events that affect cattle herds (nodes) with age structure:

Exit

Deaths or removal of cattle from a herd (n = 182,535). These events remove cattle from susceptible or infected compartments across age categories.

Enter

Births or introduction of cattle to a herd (n = 182,685). These events add susceptible cattle, typically to the youngest age category.

Internal transfer

Age transitions or within-herd movements (n = 317,081). These events move cattle between age categories within a herd, reflecting maturation and changing infection risk with age.

External transfer

Movement of cattle between herds (n = 101,472). These events transfer cattle from one herd to another across age categories, potentially introducing infected animals.

The select column in the returned data frame is mapped to the columns of the internal select matrix as follows:

• select = 1: Targets S_1 (Susceptible, age 1).

• select = 2: Targets S_2 (Susceptible, age 2).

• select = 3: Targets S_3 (Susceptible, age 3).

• select = 4: Targets S_1 and I_1 (Susceptible and Infected, age 1).

• select = 5: Targets S_2 and I_2 (Susceptible and Infected, age 2).

• select = 6: Targets S_3 and I_3 (Susceptible and Infected, age 3).

The shift column is used for Internal transfer events to define the destination compartment. It corresponds to the column index in the internal N matrix that specifies the transition (e.g., moving from age 1 to age 2).

Events are distributed across all 1,600 herds over the 4-year period. These are synthetic data generated to illustrate how to incorporate scheduled events (including births, deaths, movements, and age transitions) into an age-structured compartment model in the SimInf framework. The higher event count compared to non-age-structured models reflects the addition of internal transfer events required for age category transitions.

The data contains:

event

Event type: "exit", "enter", "intTrans", or "extTrans".

time

Day when event occurs (1-1460).

node

Affected herd identifier (1-1600).

dest

Destination herd for external transfer events, else 0.

n

Number of cattle affected.

select

Model compartment to affect (see SimInf_events).

proportion

0. Not used in this example.

shift

Determines how individuals in internal transfer events are shifted to enter another compartment.

See Also

u0_SISe3 for the corresponding initial cattle population with age structure, SISe3 for creating SISe3 models with these events and SimInf_events for event structure details

Examples

## For reproducibility, call the set.seed() function and specify the
## number of threads to use. To use all available threads, remove the
## set_num_threads() call.
set.seed(123)
set_num_threads(1)

## Create an 'SISe3' model with 1600 cattle herds (nodes) stratified
## by age, initialize it to run over 4*365 days and record data at
## weekly time-points. Add ten infected animals to age category 1 in
## the first herd to seed the outbreak.  Define 'tspan' to record the
## state of the system at weekly time-points. Load scheduled events
## events for the population of nodes with births, deaths and
## between-node movements of individuals.
u0 <- u0_SISe3
u0$I_1[1] <- 10
model <- SISe3(
    u0 = u0,
    tspan = seq(from = 1, to = 4*365, by = 7),
    events = events_SISe3,
    phi = rep(0, nrow(u0)),
    upsilon_1 = 1.8e-2,
    upsilon_2 = 1.8e-2,
    upsilon_3 = 1.8e-2,
    gamma_1 = 0.1,
    gamma_2 = 0.1,
    gamma_3 = 0.1,
    alpha = 1,
    beta_t1 = 1.0e-1,
    beta_t2 = 1.0e-1,
    beta_t3 = 1.25e-1,
    beta_t4 = 1.25e-1,
    end_t1 = 91,
    end_t2 = 182,
    end_t3 = 273,
    end_t4 = 365,
    epsilon = 0
)

## Display the number of cattle affected by each event type per day.
plot(events(model))

## Run the model to generate a single stochastic trajectory.
result <- run(model)

## Plot the median and interquartile range of the number of
## susceptible and infected individuals.
plot(result)

## Plot the proportion of nodes with at least one infected individual.
plot(result, I_1 + I_2 + I_3 ~ ., level = 2, type = "l")

## Plot the trajectory for the first herd.
plot(result, index = 1)

## Summarize the trajectory. The summary includes the number of events
## by event type.
summary(result)

Description

The global data is a numeric vector that is common to all nodes. The global data vector is passed as an argument to the transition rate functions and the post time step function.

Usage

gdata(model)

## S4 method for signature 'SimInf_model'
gdata(model)

Arguments

Value

a numeric vector

See Also

ldata for retrieving local data (node-specific), and SimInf_model for the class definition and overview of model structure.

Examples

## Create an SIR model
model <- SIR(
  u0 = data.frame(S = 99, I = 1, R = 0),
  tspan = 1:5,
  beta = 0.16,
  gamma = 0.077
)

## Set 'beta' to a new value
gdata(model, "beta") <- 2

## Extract the global data vector that is common to all nodes
gdata(model)

Description

The global data is a numeric vector that is common to all nodes. The global data vector is passed as an argument to the transition rate functions and the post time step function.

Usage

gdata(model, parameter) <- value

## S4 replacement method for signature 'SimInf_model'
gdata(model, parameter) <- value

Arguments

Value

a SimInf_model object

See Also

ldata for retrieving local data (node-specific), and SimInf_model for the class definition and overview of model structure.

Examples

## Create an SIR model
model <- SIR(
  u0 = data.frame(S = 99, I = 1, R = 0),
  tspan = 1:5,
  beta = 0.16,
  gamma = 0.077
)

## Set 'beta' to a new value
gdata(model, "beta") <- 2

## Extract the global data vector that is common to all nodes
gdata(model)

Description

Lookup individuals, in which node they are located and their age at a specified time-point.

Usage

get_individuals(x, time = NULL)

## S4 method for signature 'SimInf_individual_events'
get_individuals(x, time = NULL)

Arguments

Value

a data.frame with the columns id, node, and age.

Description

Calculate the in-degree of each node based on external transfer events ("extTrans") in the model's schedules events.

Usage

indegree(model)

Arguments

Details

The in-degree is defined as the number of unique source nodes that have sent individuals to the target node at least once. This metric measures the connectivity of the network, indicating how many different neighbors directly supply individuals to a specific node.

Value

An integer vector where each element corresponds to a node, containing the count of unique source nodes sending individuals to it.

See Also

outdegree for calculating the number of unique destination nodes a node sends individuals to. events_SIR for example event data used in network analysis. SimInf_events for details on the structure of scheduled events.

Examples

## Create an 'SIR' model with example data.
model <- SIR(
  u0 = u0_SIR(),
  tspan = 1:1460,
  events = events_SIR(),
  beta = 0.16,
  gamma = 0.077
)

## Calculate the in-degree for each node.
deg <- indegree(model)

## View the in-degree for the first 6 nodes.
head(deg)

## Plot the distribution of in-degrees across all nodes.  This
## shows how many source nodes typically send individuals to a given
## node.
hist(
  deg,
  main = "Distribution of In-Degree (Unique Sources)",
  xlab = "Number of Unique Source Nodes"
)

Description

In many countries, individual-based livestock data are collected to enable contact tracing during disease outbreaks. However, the livestock databases are not always structured in such a way that relevant information for disease spread simulations is easily retrieved. The aim of this function is to facilitate cleaning livestock event data and prepare it for usage in SimInf.

Usage

individual_events(events)

Arguments

Details

The argument events in individual_events must be a data.frame with the following columns:

• id: an integer or character identifier of the individual.

• event: four event types are supported: exit, enter, internal transfer, and external transfer. When assigning the events, they can either be coded as a numerical value or a character string: exit; 0 or 'exit', enter; 1 or 'enter', internal transfer; 2 or 'intTrans', and external transfer; 3 or 'extTrans'.

• time: an integer, character, or date (of class Date) for when the event occured. If it's a character it must be able to coerce to Date.

• node: an integer or character identifier of the source node.

• dest: an integer or character identifier of the destination node for movement events, and 'dest' will be replaced with NA for non-movement events.

Value

SimInf_individual_events

References

I. R. Ewerlöf, J. Frössling, M. Tråvén, S. Gunnarsson, L. Stengärde, E. Hurri, and S. Widgren. Exploring structural changes in the Swedish cattle population and between-holding movements. Preventive Veterinary Medicine, 243, 106608, 2025. doi:10.1016/j.prevetmed.2025.106608

See Also

node_events.

Description

W0(x) is the principal branch of the solution of the function defined by $We^W = x$ for $x >= -1/e$. The value is calculated using GNU Scientific Library (GSL).

Usage

lambertW0(x)

Arguments

References

GNU Scientific Library <https://www.gnu.org/software/gsl/>

Examples

## Should equal 1, as 1 * exp(1) = e.
lambertW0(exp(1))

## Should equal 0, as 0 * exp(0) = 0.
lambertW0(0)

## Should equal -1.
lambertW0(-exp(-1))

Description

The local data is a numeric vector that is specific to a node. The local data vector is passed as an argument to the transition rate functions and the post time step function.

Usage

ldata(model, node)

## S4 method for signature 'SimInf_model'
ldata(model, node)

Arguments

Value

a numeric vector

See Also

gdata for retrieving global data (common to all nodes), and SimInf_model for the class definition and overview of model structure.

Examples

## Create an 'SISe' model with 1600 nodes.
model <- SISe(
  u0 = u0_SIS(),
  tspan = 1:100,
  events = events_SIS(),
  phi = 0,
  upsilon = 1.8e-2,
  gamma = 0.1,
  alpha = 1,
  beta_t1 = 1.0e-1,
  beta_t2 = 1.0e-1,
  beta_t3 = 1.25e-1,
  beta_t4 = 1.25e-1,
  end_t1 = c(91, 101),
  end_t2 = c(182, 185),
  end_t3 = c(273, 275),
  end_t4 = c(365, 360),
  epsilon = 0
)

## Display local data from the first two nodes.
ldata(model, node = 1)
ldata(model, node = 2)

Description

Get the number of iterations (samples) in the Markov Chain Monte Carlo (MCMC) chain stored in a SimInf_pmcmc object.

Usage

## S4 method for signature 'SimInf_pmcmc'
length(x)

Arguments

Value

An integer scalar representing the number of rows in the chain slot (i.e., the total number of samples in the chain).

Description

Extract the estimated log likelihood from a SimInf_pfilter object.

Usage

## S4 method for signature 'SimInf_pfilter'
logLik(object)

Arguments

Value

the estimated log likelihood.

Description

Describe your model using a simple string syntax in R. mparse parses this description, generates model-specific C code, and returns a SimInf_model object ready for simulation.

Usage

mparse(
  transitions = NULL,
  compartments = NULL,
  ldata = NULL,
  gdata = NULL,
  u0 = NULL,
  v0 = NULL,
  tspan = NULL,
  events = NULL,
  E = NULL,
  N = NULL,
  pts_fun = NULL,
  use_enum = FALSE,
  pre_code = NULL
)

Arguments

Value

a SimInf_model object

See Also

SimInf_model for the class definition of the returned model object. SIR, SEIR, SIS, SISe for high-level model constructors that use predefined structures. Vignette "Getting started with mparse" for a comprehensive tutorial on defining custom models, including syntax, variables, and event handling. package_skeleton for creating an installable R package from a mparse model.

Examples

## For reproducibility, set the seed.
set.seed(22)

## Create an SIR model with a defined population size variable.
model <- mparse(
  transitions = c(
    "S -> beta * S * I / N -> I",
    "I -> gamma * I -> R",
    "N <- S + I + R"
  ),
  compartments = c("S", "I", "R"),
  gdata = c(beta = 0.16, gamma = 0.077),
  u0 = data.frame(S = 100, I = 1, R = 0),
  tspan = 1:100
)

## Run and plot the result.
result <- run(model)
plot(result)

Description

Extract the number of compartments from a SimInf_model object. Compartments represent the distinct states an individual can occupy within a node (e.g., Susceptible, Infected, Recovered). This count is equivalent to the number of columns in the initial state vector u0.

Usage

n_compartments(model)

## S4 method for signature 'SimInf_model'
n_compartments(model)

Arguments

Value

An integer scalar representing the total number of compartments in the model.

Examples

## Create an 'SIR' model with 3 compartments (S, I, R).
u0 <- data.frame(
  S = rep(99, 100),
  I = rep(1, 100),
  R = rep(0, 100)
)

model <- SIR(
  u0 = u0,
  tspan = 1:10,
  beta = 0.16,
  gamma = 0.077
)

## Get the number of compartments.
n_compartments(model)

Description

Extract the number of generations performed in an Approximate Bayesian Computation (ABC) analysis from a SimInf_abc object. Each generation represents a sequential round of the algorithm, involving the simulation of particles, acceptance/rejection based on the distance threshold, and parameter perturbation for the next round.

Usage

n_generations(object)

## S4 method for signature 'SimInf_abc'
n_generations(object)

Arguments

Value

An integer scalar representing the total number of generations executed in the analysis.

Description

Extract the number of nodes from a SimInf_model object. A node represents a distinct sub-population in the model, but its definition is determined by the modeller. For example, a node can represent a cattle herd, a pen within a herd, or even a single individual, depending on the research question and the scale of the study. This count is equivalent to the number of rows in the initial state vector u0.

Usage

n_nodes(model)

## S4 method for signature 'SimInf_model'
n_nodes(model)

## S4 method for signature 'SimInf_pfilter'
n_nodes(model)

## S4 method for signature 'SimInf_pmcmc'
n_nodes(model)

Arguments

Value

An integer scalar representing the total number of nodes in the model.

Examples

## Create an 'SIR' model with 100 nodes.
u0 <- data.frame(
  S = rep(99, 100),
  I = rep(1, 100),
  R = rep(0, 100)
)

model <- SIR(
  u0 = u0,
  tspan = 1:10,
  beta = 0.16,
  gamma = 0.077
)

## Get the number of nodes.
n_nodes(model)

Description

Extract the number of replicates from a SimInf_model object. Replicates are independent copies of the model state used by the particle filter (pfilter). This value is set internally by pfilter and is not specified when creating a standard model. If the model has not been processed by a filter, the value will be 1.

Usage

n_replicates(model)

## S4 method for signature 'SimInf_model'
n_replicates(model)

Arguments

Value

An integer scalar representing the number of replicates in the model.

Examples

## Create a standard 'SIR' model.
u0 <- data.frame(
  S = rep(99, 100),
  I = rep(1, 100),
  R = rep(0, 100)
)

model <- SIR(
  u0 = u0,
  tspan = 1:10,
  beta = 0.16,
  gamma = 0.077
)

## Get the number of replicates (default is 1 for standard
## models).
n_replicates(model)

Description

In many countries, individual-based livestock data are collected to enable contact tracing during disease outbreaks. However, the livestock databases are not always structured in such a way that relevant information for disease spread simulations is easily retrieved. The aim of this function is to facilitate cleaning livestock event data and prepare it for usage in SimInf.

Usage

node_events(x, time = NULL, target = NULL, age = NULL)

## S4 method for signature 'SimInf_individual_events'
node_events(x, time = NULL, target = NULL, age = NULL)

Arguments

Details

The individual-based events will be aggregated on node-level. The select value is determined by the event type and age category. If there is only one age category, i.e., age=NULL, then select=1 for the enter events, and select=2 for all other events. If there are two age categories, then select=1 for the enter events in the first age category, and select=2 for the enter events in the second age category. Similarly, select=3 for all other events in the first age category, and select=4 for all other events in the first second category. With three age categories, it works similarly with select=1,2,3 for the enter events in each age category, respectively. And select=4,5,6 for all other events.

Value

a data.frame with the columns event, time, node, dest, n, proportion, select, and shift.

See Also

individual_events.

Description

Example data to initialize a population of 1600 nodes and demonstrate various models.

Usage

data(nodes)

Format

A data.frame

Examples

## Not run: 
## For reproducibility, call the set.seed() function and specify
## the number of threads to use. To use all available threads,
## remove the set_num_threads() call.
set.seed(123)
set_num_threads(1)

## Create an 'SIR' model with 1600 nodes and initialize
## it to run over 4*365 days. Add one infected individual
## to the first node.
u0 <- u0_SIR()
u0$I[1] <- 1
tspan <- seq(from = 1, to = 4*365, by = 1)
model <- SIR(u0     = u0,
             tspan  = tspan,
             events = events_SIR(),
             beta   = 0.16,
             gamma  = 0.077)

## Run the model to generate a single stochastic trajectory.
result <- run(model)

## Determine nodes with one or more infected individuals in the
## trajectory. Extract the 'I' compartment and check for any
## infected individuals in each node.
infected <- colSums(trajectory(result, ~ I, format = "matrix")) > 0

## Display infected nodes in 'blue' and non-infected nodes in 'yellow'.
data("nodes", package = "SimInf")
col <- ifelse(infected, "blue", "yellow")
plot(y ~ x, nodes, col = col, pch = 20, cex = 2)

## End(Not run)

Description

Calculate the out-degree of each node based on external transfer events ("extTrans") in the model's scheduled events.

Usage

outdegree(model)

Arguments

Details

The out-degree is defined as the number of unique destination nodes that receive individuals from the source node at least once. This metric measures the connectivity of the network, indicating how many different neighbors a specific node directly sends individuals to.

Value

An integer vector where each element corresponds to a node, containing the count of unique destination nodes receiving individuals from it.

See Also

indegree for calculating the number of unique source nodes that send individuals to a node. events_SIR for example event data used in network analysis. SimInf_events for details on the structure of scheduled events.

Examples

## Create an 'SIR' model with example data.
model <- SIR(
  u0 = u0_SIR(),
  tspan = 1:1460,
  events = events_SIR(),
  beta = 0.16,
  gamma = 0.077
)

## Calculate the out-degree for each node.
deg <- outdegree(model)

## View the out-degree for the first 6 nodes.
head(deg)

## Plot the distribution of out-degrees across all nodes.
## This shows how many destination nodes typically receive
## individuals from a given node.
hist(
  deg,
  main = "Distribution of Out-Degree (Unique Destinations)",
  xlab = "Number of Unique Destination Nodes"
)

Description

Generate a source package directory from a SimInf_model object, allowing the model to be installed as a standalone add-on R package. This is useful for distributing custom models or for improving startup performance by compiling the C code once during package installation rather than on the first call to run().

Usage

package_skeleton(
  model,
  name = NULL,
  path = ".",
  author = NULL,
  email = NULL,
  maintainer = NULL,
  license = "GPL-3"
)

Arguments

Details

The typical workflow is:

1. Define the model using mparse or SimInf_model.

2. Call package_skeleton to generate the package

3. Install the package using install.packages or R CMD INSTALL.

Value

invisible(NULL).

References

Read the Writing R Extensions manual for more details on building R packages.

See Also

mparse for creating a SimInf_model object from a model specification. INSTALL and install.packages for installing the generated package.

Description

Produce a matrix of scatterplots showing the relationship between the number of individuals in different model compartments. The ijth panel displays the counts of compartment j plotted against compartment i.

Usage

## S4 method for signature 'SimInf_model'
pairs(x, compartments = NULL, index = NULL, ...)

Arguments

Details

Data is aggregated across all time points in tspan and the selected nodes (specified by index). This allows for visualizing the correlation structure between compartments throughout the simulation.

Examples

## For reproducibility, set the seed and number of threads.
set.seed(123)
set_num_threads(1)

## Create an 'SIR' model with 10 nodes.
model <- SIR(
  u0 = data.frame(
    S = rep(99, 10),
    I = rep(1, 10),
    R = rep(0, 10)
  ),
  tspan = 1:100,
  beta = 0.16,
  gamma = 0.077
)

## Run the model.
result <- run(model)

## Create a scatterplot matrix for all compartments across all
## nodes.
pairs(result)

## Create a scatterplot matrix for the S and I compartments,
## using only data from nodes 1 and 2.
pairs(result, compartments = c("S", "I"), index = 1:2)

Description

The bootstrap filtering algorithm. Systematic resampling is performed at each observation.

Usage

pfilter(model, obs_process, data, n_particles, init_model = NULL)

## S4 method for signature 'SimInf_model'
pfilter(model, obs_process, data, n_particles, init_model = NULL)

Arguments

Value

A SimInf_pfilter object.

References

N. J. Gordon, D. J. Salmond, and A. F. M. Smith. Novel Approach to Nonlinear/Non-Gaussian Bayesian State Estimation. Radar and Signal Processing, IEE Proceedings F, 140(2) 107–113, 1993. doi:10.1049/ip-f-2.1993.0015

Examples

## Not run: 
## Let us consider an SIR model in a closed population with N = 100
## individuals of whom one is initially infectious and the rest are
## susceptible. First, generate one realisation (with a specified
## seed) from the model with known parameters 'beta = 0.16' and
## 'gamma = 0.077'. Then, use 'pfilter' to apply the bootstrap
## particle algorithm on the simulated data.
model <- SIR(u0 = data.frame(S = 99, I = 1, R = 0),
             tspan = seq(1, 100, by = 3),
             beta = 0.16,
             gamma = 0.077)

## Run the SIR model to generate simulated observed data for the
## number of infected individuals.
set.seed(22)
infected <- trajectory(run(model), "I")[, c("time", "I")]
colnames(infected) <- c("time", "Iobs")

## Use a Poison observation process for the infected individuals, such
## that 'Iobs ~ poison(I + 1e-6)'. A small constant '1e-6' is added to
## prevent numerical errors, since the simulated counts 'I' could be
## zero, which would result in the Poisson rate parameter being zero,
## which violates the conditions of the Poisson distribution. Use 1000
## particles.
pf <- pfilter(model,
              obs_process = Iobs ~ poisson(I + 1e-6),
              data = infected,
              n_particles = 1000)

## Print a brief summary.
pf

## Compare the number infected 'I' in the filtered trajectory with the
## infected 'Iobs' in the observed data.
plot(pf, ~I)
lines(Iobs ~ time, infected, col = "blue", lwd = 2, type = "s")

## End(Not run)

Description

Produce diagnostic plots of the Approximate Bayesian Computation (ABC) posterior distribution stored in a SimInf_abc object.

Usage

## S4 method for signature 'SimInf_abc'
plot(x, y, ...)

Arguments

Details

The function generates a scatterplot matrix of the parameter values for the specified generation:

• Diagonal panels: Display a normalized density estimate with a rug plot showing individual samples.

• Upper triangular panels: Display scatterplots of the raw parameter samples for each pair of variables.

• Lower triangular panels: Display contour lines representing the 2D kernel density estimate of the joint distribution between parameter pairs.

If only a single parameter is selected, a single density plot with a rug is produced.

Description

Display the distribution of scheduled events over time

Usage

## S4 method for signature 'SimInf_events'
plot(x, frame.plot = FALSE, ...)

Arguments

Description

Display the distribution of individual events over time

Usage

## S4 method for signature 'SimInf_individual_events'
plot(x, frame.plot = FALSE, ...)

Arguments

Description

Plot the median and quantile range of the counts in all nodes, the counts in specified nodes, or the prevalence of a disease. The function supports formula notation for specifying compartments and prevalence calculations.

Usage

## S4 method for signature 'SimInf_model'
plot(
  x,
  y,
  level = 1,
  index = NULL,
  range = 0.5,
  type = "s",
  lwd = 2,
  frame.plot = FALSE,
  legend = TRUE,
  log = "",
  ...
)

Arguments

Details

For a comprehensive tutorial with detailed explanations and more examples, see the vignette("Post-process data in a trajectory", package = "SimInf").

Examples

## For reproducibility, set the seed and number of threads.
set.seed(123)
set_num_threads(1)

## Create an 'SIR' model with 100 nodes.
model <- SIR(u0 = data.frame(S = rep(990, 100),
                             I = rep(10, 100),
                             R = rep(0, 100)),
             tspan = 1:100,
             beta = 0.16,
             gamma = 0.077)

## Run the model.
result <- run(model)

## Plot counts (median and IQR)
plot(result)

## Plot individual trajectories for specific nodes
plot(result, index = 1:3, range = FALSE)

## Plot prevalence (proportion of infected)
plot(result, I ~ S + I + R)

## Customize labels
plot(result, "I", xlab = "Time", ylab = "Count", main = "Infections")

Description

Diagnostic plot of a particle filter object

Usage

## S4 method for signature 'SimInf_pfilter'
plot(x, y, ...)

Arguments

Description

Display the (approximate) posterior distributions obtained from fitting a particle Markov chain Monte Carlo algorithm, or the corresponding trace plots.

Usage

## S4 method for signature 'SimInf_pmcmc'
plot(x, y, start = 1, end = NULL, thin = 1, ...)

Arguments

Description

Estimates model parameters using the PMCMC algorithm, which combines a particle filter for likelihood evaluation with a Metropolis-Hastings MCMC sampler for parameter estimation.

Usage

pmcmc(
  model,
  obs_process,
  data,
  priors,
  n_particles,
  n_iterations,
  theta = NULL,
  covmat = NULL,
  adaptmix = 0.05,
  adaptive = 100,
  post_proposal = NULL,
  init_model = NULL,
  post_particle = NULL,
  chain = NULL,
  verbose = FALSE
)

## S4 method for signature 'SimInf_model'
pmcmc(
  model,
  obs_process,
  data,
  priors,
  n_particles,
  n_iterations,
  theta = NULL,
  covmat = NULL,
  adaptmix = 0.05,
  adaptive = 100,
  post_proposal = NULL,
  init_model = NULL,
  post_particle = NULL,
  chain = NULL,
  verbose = FALSE
)

Arguments

Value

A SimInf_pmcmc object containing the fitted parameters and diagnostic information for all iterations.

References

C. Andrieu, A. Doucet and R. Holenstein. Particle Markov chain Monte Carlo methods. Journal of the Royal Statistical Society, Series B 72, 269–342, 2010. doi:10.1111/j.1467-9868.2009.00736.x

G. O. Roberts and J. S. Rosenthal. Examples of adaptive MCMC. Journal of computational and graphical statistics, 18(2), 349–367, 2009. doi:10.1198/jcgs.2009.06134

See Also

continue_pmcmc for running additional iterations.

Description

Calculate the proportion of cases (specified on the left-hand side of the formula) relative to the population at risk (specified on the right-hand side) at different aggregation levels. The function supports three levels of calculation:

• Level 1 (Population Prevalence): The proportion of cases in the total population at risk across all nodes.

• Level 2 (Node Prevalence): The proportion of nodes that contain at least one case, considering only nodes where the population at risk is greater than zero.

• Level 3 (Within-Node Prevalence): The proportion of cases relative to the population at risk calculated separately for each node.

Usage

prevalence(model, formula, level = 1, index = NULL, ...)

Arguments

Description

Calculate the proportion of cases (specified on the left-hand side of the formula) relative to the population at risk (specified on the right-hand side) at different aggregation levels. The function supports three levels of calculation:

• Level 1 (Population Prevalence): The proportion of cases in the total population at risk across all nodes.

• Level 2 (Node Prevalence): The proportion of nodes that contain at least one case, considering only nodes where the population at risk is greater than zero.

• Level 3 (Within-Node Prevalence): The proportion of cases relative to the population at risk calculated separately for each node.

Usage

## S4 method for signature 'SimInf_model'
prevalence(model, formula, level, index, format = c("data.frame", "matrix"))

Arguments

Value

A data.frame if format = "data.frame", else a matrix.

Examples

## For reproducibility, set the seed and number of threads.
set.seed(1)
set_num_threads(1)

## Create an 'SIR' model with 6 nodes.
u0 <- data.frame(
  S = 100:105,
  I = c(0, 1, 0, 2, 0, 3),
  R = rep(0, 6)
)

model <- SIR(
  u0 = u0,
  tspan = 1:10,
  beta = 0.16,
  gamma = 0.077
)

## Run the model.
result <- run(model)

## 1. Population Prevalence (level = 1, default)
## Proportion of infected individuals in the total population.
prevalence(result, I ~ S + I + R)

## Shorthand: '.' represents all compartments in the model (S + I + R).
prevalence(result, I ~ .)

## 2. Node Prevalence (level = 2)
## Proportion of nodes with at least one infected individual.
prevalence(result, I ~ S + I + R, level = 2)

## 3. Within-Node Prevalence (level = 3)
## Prevalence calculated separately for each node.
prevalence(result, I ~ S + I + R, level = 3)

## 4. Conditional Prevalence
## Calculate prevalence only in nodes where the number of
## recovered is zero.
prevalence(result, I ~ S + I + R | R == 0, level = 3)

Description

Calculate the proportion of cases (specified on the left-hand side of the formula) relative to the population at risk (specified on the right-hand side) at different aggregation levels. The function supports three levels of calculation:

• Level 1 (Population Prevalence): The proportion of cases in the total population at risk across all nodes.

• Level 2 (Node Prevalence): The proportion of nodes that contain at least one case, considering only nodes where the population at risk is greater than zero.

• Level 3 (Within-Node Prevalence): The proportion of cases relative to the population at risk calculated separately for each node.

Usage

## S4 method for signature 'SimInf_pfilter'
prevalence(model, formula, level, index, format = c("data.frame", "matrix"))

Arguments

Value

A data.frame if format = "data.frame", else a matrix.

Description

Calculate the proportion of cases (specified on the left-hand side of the formula) relative to the population at risk (specified on the right-hand side) at different aggregation levels. The function supports three levels of calculation:

• Level 1 (Population Prevalence): The proportion of cases in the total population at risk across all nodes.

• Level 2 (Node Prevalence): The proportion of nodes that contain at least one case, considering only nodes where the population at risk is greater than zero.

• Level 3 (Within-Node Prevalence): The proportion of cases relative to the population at risk calculated separately for each node.

Usage

## S4 method for signature 'SimInf_pmcmc'
prevalence(model, formula, level, index, start = 1, end = NULL, thin = 1)

Arguments

Value

A data.frame where the first column is the iteration and the remaining columns are the result from calling prevalence with the arguments formula, level and index for each iteration.

Description

Define a template (or "punchcard") specifying which data points (nodes, time-points, and compartments) should be recorded during a simulation. This feature is useful for saving memory when the model has many nodes or time-points, but only a subset of the data is needed for post-processing.

Usage

punchcard(model) <- value

## S4 replacement method for signature 'SimInf_model'
punchcard(model) <- value

Arguments

Details

The template is specified as a data.frame with columns:

• time: The time-points to record.

• node: The node indices to record.

• CompartmentName: Logical columns (e.g., S, I, R) indicating whether to record that specific compartment for the corresponding row. If a compartment column is omitted, it defaults to FALSE (not recorded). If only time and node are provided, all compartments are recorded for those points.

Important: The specified time values, node indices, and compartment names must exist in the model. If any value in the template does not match the model's configuration (e.g., a time point not in tspan or a compartment not defined in the model), an error will be raised when the template is set.

Note that the template only affects which data is stored in the result object; it does not affect how the solver simulates the trajectory.

Examples

## For reproducibility, set the seed and number of threads.
set.seed(123)
set_num_threads(1)

## Create an 'SIR' model with 6 nodes
u0 <- data.frame(
  S = 100:105,
  I = 1:6,
  R = rep(0, 6)
)
model <- SIR(
  u0 = u0,
  tspan = 1:10,
  beta = 0.16,
  gamma = 0.077
)

## Run the model with default recording (all data)
result_full <- run(model)
head(trajectory(result_full))

## Define a template to record only nodes 2 and 4 at times 3 and 5
## for all compartments.
df <- data.frame(
  time = c(3, 5, 3, 5),
  node = c(2, 2, 4, 4),
  S = TRUE, I = TRUE, R = TRUE
)
punchcard(model) <- df

result_sparse <- run(model)
trajectory(result_sparse)

## Record only specific compartments (e.g., S and R, but not I)
df <- data.frame(
  time = c(3, 5, 3, 5),
  node = c(2, 2, 4, 4),
  S = TRUE, I = FALSE, R = TRUE
)
punchcard(model) <- df
result_partial <- run(model)
trajectory(result_partial)

## Shortcut: If only 'time' and 'node' are specified, all
## compartments are recorded.
df <- data.frame(
  time = c(3, 5, 3, 5),
  node = c(2, 2, 4, 4)
)
punchcard(model) <- df

## Reset to record all data (equivalent to no template)
punchcard(model) <- NULL
result_reset <- run(model)
head(trajectory(result_reset))

Description

Run the SimInf stochastic simulation algorithm

Usage

run(model, ...)

## S4 method for signature 'SimInf_model'
run(model, solver = c("ssm", "aem"), ...)

## S4 method for signature 'SEIR'
run(model, solver = c("ssm", "aem"), ...)

## S4 method for signature 'SIR'
run(model, solver = c("ssm", "aem"), ...)

## S4 method for signature 'SIS'
run(model, solver = c("ssm", "aem"), ...)

## S4 method for signature 'SISe'
run(model, solver = c("ssm", "aem"), ...)

## S4 method for signature 'SISe3'
run(model, solver = c("ssm", "aem"), ...)

## S4 method for signature 'SISe3_sp'
run(model, solver = c("ssm", "aem"), ...)

## S4 method for signature 'SISe_sp'
run(model, solver = c("ssm", "aem"), ...)

## S4 method for signature 'SimInf_abc'
run(model, ...)

Arguments

Value

SimInf_model object with result from simulation.

References

S. Widgren, P. Bauer, R. Eriksson and S. Engblom. SimInf: An R Package for Data-Driven Stochastic Disease Spread Simulations. Journal of Statistical Software, 91(12), 1–42, 2019. doi:10.18637/jss.v091.i12. An updated version of this paper is available as a vignette in the package.

P. Bauer, S. Engblom and S. Widgren. Fast Event-Based Epidemiological Simulations on National Scales. International Journal of High Performance Computing Applications, 30(4), 438–453, 2016. doi: 10.1177/1094342016635723

P. Bauer and S. Engblom. Sensitivity Estimation and Inverse Problems in Spatial Stochastic Models of Chemical Kinetics. In: A. Abdulle, S. Deparis, D. Kressner, F. Nobile and M. Picasso (eds.), Numerical Mathematics and Advanced Applications - ENUMATH 2013, pp. 519–527, Lecture Notes in Computational Science and Engineering, vol 103. Springer, Cham, 2015. doi:10.1007/978-3-319-10705-9_51

Examples

## For reproducibility, call the set.seed() function and specify
## the number of threads to use. To use all available threads,
## remove the set_num_threads() call.
set.seed(123)
set_num_threads(1)

## Create an 'SIR' model with 10 nodes and initialise
## it to run over 100 days.
model <- SIR(u0 = data.frame(S = rep(99, 10),
                             I = rep(1, 10),
                             R = rep(0, 10)),
             tspan = 1:100,
             beta = 0.16,
             gamma = 0.077)

## Run the model and save the result.
result <- run(model)

## Plot the proportion of susceptible, infected and recovered
## individuals.
plot(result)

Description

Create an SEIR model to be used by the simulation framework.

Usage

SEIR(u0, tspan, events = NULL, beta = NULL, epsilon = NULL, gamma = NULL)

Arguments

Details

The SEIR model extends the standard SIR model by adding an Exposed (E) compartment for individuals who have been infected but are not yet infectious. This accounts for the latent period of the disease.

The model is defined by three state transitions:

$S \stackrel{\beta S I / N}{\longrightarrow} E$

$E \stackrel{\epsilon E}{\longrightarrow} I$

$I \stackrel{\gamma I}{\longrightarrow} R$

where $\beta$ is the transmission rate, $\epsilon$ is the incubation rate (inverse of the latent period), $\gamma$ is the recovery rate, and $N = S + E + I + R$ is the total population size in each node. Here, $S$, $E$, $I$, and $R$ represent the number of susceptible, exposed, infected, and recovered individuals in that specific node.

The argument u0 must be a data.frame with one row for each node with the following columns:

S

The number of susceptible individuals in each node

E

The number of exposed individuals in each node

I

The number of infected individuals in each node

R

The number of recovered individuals in each node

Value

A SimInf_model of class SEIR

See Also

SEIR for the class definition. SIR, SIS, SISe, SISe3 and SISe_sp for other predefined models. mparse for creating custom models. run for running the simulation. trajectory, prevalence and plot for post-processing and visualization.

Examples

## For reproducibility, set the seed.
set.seed(3)

## Create an SEIR model object.
model <- SEIR(
  u0 = data.frame(S = 99, E = 0, I = 1, R = 0),
  tspan = 1:100,
  beta = 0.16,
  epsilon = 0.25,
  gamma = 0.077
)

## Run the SEIR model and plot the result.
result <- run(model)
plot(result)

Description

Class to handle the SEIR model. This class inherits from SimInf_model, meaning that SEIR objects are fully compatible with all generic functions defined for SimInf_model, such as run, plot, trajectory, and prevalence.

Details

The SEIR model extends the standard SIR model by adding an Exposed (E) compartment for individuals who have been infected but are not yet infectious. This accounts for the latent period of the disease.

The model is defined by three state transitions:

$S \stackrel{\beta S I / N}{\longrightarrow} E$

$E \stackrel{\epsilon E}{\longrightarrow} I$

$I \stackrel{\gamma I}{\longrightarrow} R$

where $\beta$ is the transmission rate, $\epsilon$ is the incubation rate (inverse of the latent period), $\gamma$ is the recovery rate, and $N = S + E + I + R$ is the total population size in each node. Here, $S$, $E$, $I$, and $R$ represent the number of susceptible, exposed, infected, and recovered individuals in that specific node.

See Also

SEIR for creating an SEIR model object, SimInf_model for the parent class definition, and SIR for the base model without the latent period.

Description

Utility function to extract events@E from a SimInf_model object, see SimInf_events

Usage

select_matrix(model)

## S4 method for signature 'SimInf_model'
select_matrix(model)

Arguments

Value

dgCMatrix object.

Examples

## Create an SIR model
model <- SIR(u0 = data.frame(S = 99, I = 1, R = 0),
             tspan = 1:5, beta = 0.16, gamma = 0.077)

## Extract the select matrix from the model
select_matrix(model)

Description

Utility function to set events@E in a SimInf_model object, see SimInf_events.

Usage

select_matrix(model) <- value

## S4 replacement method for signature 'SimInf_model'
select_matrix(model) <- value

Arguments

Examples

## Create an SIR model.
model <- SIR(u0 = data.frame(S = 99, I = 1, R = 0),
             tspan = 1:5, beta = 0.16, gamma = 0.077)

## Set the select matrix.
select_matrix(model) <- matrix(c(1, 0, 0, 1, 1, 1, 0, 0, 1), nrow = 3)

## Extract the select matrix from the model.
select_matrix(model)

## Set the select matrix using a data.frame instead.
select_matrix(model) <- data.frame(
    compartment = c("S", "S", "I", "R", "R"),
    select      = c( 1,   2,   2,   2,   3))

## Extract the select matrix from the model.
select_matrix(model)

Description

Set the number of threads to be used in SimInf code that is parallelized with OpenMP (if available).

Usage

set_num_threads(threads = NULL)

Arguments

Details

The number of threads is initialized when SimInf is first loaded in the R session, based on optional environment variables (see ‘Details’). It can also be explicitly set by calling set_num_threads. If the environment variables affecting the thread count change, set_num_threads must be called again for the new values to take effect.

The function determines the number of available processors using omp_get_num_procs() and limits the thread count based on omp_get_thread_limit() and the following environment variables (checked in order of precedence):

• SIMINF_NUM_THREADS: Specific limit for SimInf.

• OMP_NUM_THREADS: General OpenMP limit.

• OMP_THREAD_LIMIT: Maximum thread limit.

The threads argument allows you to override these limits, provided the requested value does not exceed the maximum allowed by the environment or hardware constraints.

Value

The previous value of the thread count is returned invisibly.

Description

Utility function to extract the shift matrix events@N from a SimInf_model object, see SimInf_events

Usage

shift_matrix(model)

## S4 method for signature 'SimInf_model'
shift_matrix(model)

Arguments

Value

A mtrix.

Examples

## Create an SIR model
model <- SIR(u0 = data.frame(S = 99, I = 1, R = 0),
             tspan = 1:5, beta = 0.16, gamma = 0.077)

## Extract the shift matrix from the model
shift_matrix(model)

Description

Utility function to set events@N in a SimInf_model object, see SimInf_events.

Usage

shift_matrix(model) <- value

## S4 replacement method for signature 'SimInf_model'
shift_matrix(model) <- value

Arguments

Value

SimInf_model object

Examples

## Create an SIR model
model <- SIR(u0 = data.frame(S = 99, I = 1, R = 0),
             tspan = 1:5, beta = 0.16, gamma = 0.077)

## Set the shift matrix.
shift_matrix(model) <- matrix(c(2, 1, 0), nrow = 3)

## Extract the shift matrix from the model.
shift_matrix(model)

## Set the shift matrix using a data.frame instead.
shift_matrix(model) <- data.frame(
    compartment = c("S", "I", "R"),
    shift = c(1, 1, 1),
    value = c(2, 1, 0))

## Extract the shift matrix from the model.
shift_matrix(model)

Description

Print summary of a SimInf_abc object

Usage

## S4 method for signature 'SimInf_abc'
show(object)

Arguments

Value

invisible(object).

Description

Shows the number of scheduled events.

Usage

## S4 method for signature 'SimInf_events'
show(object)

Arguments

Value

None (invisible 'NULL').

Description

Print summary of a SimInf_individual_events object

Usage

## S4 method for signature 'SimInf_individual_events'
show(object)

Arguments

Value

invisible(object).

Description

Brief summary of SimInf_model

Usage

## S4 method for signature 'SimInf_model'
show(object)

Arguments

Value

None (invisible 'NULL').

Examples

## Create an 'SIR' model with 10 nodes and initialise
## it to run over 100 days.
model <- SIR(u0 = data.frame(S = rep(99, 10),
                             I = rep(1, 10),
                             R = rep(0, 10)),
             tspan = 1:100,
             beta = 0.16,
             gamma = 0.077)

## Brief summary of the model
model

## Run the model and save the result
result <- run(model)

## Brief summary of the result. Note that 'U' and 'V' are
## non-empty after running the model.
result

Description

Brief summary of a SimInf_pfilter object

Usage

## S4 method for signature 'SimInf_pfilter'
show(object)

Arguments

Value

invisible(object).

Description

Brief summary of a SimInf_pmcmc object

Usage

## S4 method for signature 'SimInf_pmcmc'
show(object)

Arguments

Value

invisible(object).

Description

The SimInf package provides a flexible framework for data-driven spatio-temporal disease spread modeling, combining the speed of compiled C code with the flexibility of R. It is designed to efficiently simulate disease transmission dynamics alongside population demographics and dynamic contact networks.

Details

The SimInf framework models infection dynamics within each subpopulation (node) as continuous-time Markov chains (CTMC) using the Gillespie stochastic simulation algorithm (SSA). Additionally, SimInf can incorporate data—such as births, deaths, and movements—as scheduled events. These events trigger at predefined time points and modify the state of subpopulations by randomly selecting a specified number of individuals from the affected compartments. This capability allows simulations to be driven by empirical records or synthetic scenarios while maintaining the stochastic nature of population demographics and movement networks.

The package supports both predefined models (e.g., SIR, SIS) and custom model specifications via the mparse function. mparse serves as the primary interface for defining custom compartment models, allowing users to describe transitions using a simple, human-readable string syntax in R. The function then parses this description, generates model-specific C code, and returns a SimInf_model object ready for simulation. This approach combines the flexibility of R with the computational speed of compiled code, making it well-suited for models with complex propensity functions, multiple compartments, or node-specific parameters. See the vignette "Getting started with mparse" for a detailed tutorial on defining custom models.

To facilitate the use of real livestock data, SimInf provides functionality for preparing individual event records (e.g., births, deaths, and movements) for simulation. The individual_events function cleans and validates raw livestock events, resolving inconsistencies and structuring the data into the format required by the simulation. This process transforms complex, individual-level logs into the scheduled events that drive population demographics and movement networks. See the vignette "Scheduled events" for a detailed guide on processing livestock data and integrating it into disease spread models.

After a model is created, a simulation is executed using the run method. Upon successful completion, run returns a new SimInf_model object containing the original configuration plus the simulated stochastic trajectory.

To inspect and analyze the results, SimInf provides a suite of utility functions:

• summary and show for a quick overview of the model structure and simulation results.

• plot for visualizing the time series of compartments and continuous state variables.

• trajectory for extracting the full time series data for custom analysis.

• prevalence for calculating and summarizing disease prevalence across nodes and time.

These functions facilitate both rapid exploratory analysis and detailed post-processing of simulation outcomes. See the vignette "Post-process data in a trajectory" for a comprehensive tutorial on extracting and analyzing simulation results.

Beyond simulation, the package provides functionality to fit models to time series data using two Bayesian inference methods:

• Approximate Bayesian Computation Sequential Monte Carlo (ABC-SMC), implemented in abc, based on the approach by Toni and others (2009) doi:10.1098/rsif.2008.0172.

• Particle Markov Chain Monte Carlo (PMCMC), implemented in pmcmc, based on the approach by Andrieu and others (2010) doi:10.1111/j.1467-9868.2009.00736.x.

Both methods enable parameter estimation in stochastic models where the likelihood function is intractable, by using simulated data to estimate the posterior distributions of model parameters.

Finally, one of the core design goals of SimInf is extensibility. The package provides functionality to generate the required C and R code from a model specification, enabling users to create custom R packages that leverage the numerical solvers of SimInf. This facilitates complex epidemiological research by allowing researchers to encapsulate their specific model structures within a standard R package, ensuring reproducibility and ease of sharing. See package_skeleton for details on creating a new R package based on a custom model specification.

Author(s)

Maintainer: Stefan Widgren [email protected] (ORCID)

Authors:

• Stefan Widgren [email protected] (ORCID)

• Robin Eriksson (ORCID)

• Stefan Engblom (ORCID)

• Pavol Bauer (ORCID)

Other contributors:

• Thomas Rosendal (ORCID) [contributor]

• Ivana Rodriguez Ewerlöf (ORCID) [contributor]

• Attractive Chaos (Author of 'kvec.h'.) [copyright holder]

References

S. Widgren, P. Bauer, R. Eriksson and S. Engblom. SimInf: An R Package for Data-Driven Stochastic Disease Spread Simulations. Journal of Statistical Software, 91(12), 1–42, 2019. doi:10.18637/jss.v091.i12. An updated version of this paper is available as a vignette in the package.

See Also

Useful links:

• https://github.com/stewid/SimInf

• https://stewid.github.io/SimInf/

• Report bugs at https://github.com/stewid/SimInf/issues

Description

Storage class for the results of an Approximate Bayesian Computation (ABC) parameter estimation using Sequential Monte Carlo (SMC). The SimInf_abc class holds the model definition, prior distributions, accepted parameter values (particles), weights, distances, and convergence diagnostics.

Slots

model

A SimInf_model object containing the model structure (transitions, compartments, etc.) for which parameters are being estimated.

priors

A data.frame defining the prior distributions for the parameters. It contains four columns:

• parameter: The name of the parameter in the model.

• distribution: The prior distribution type. Valid values are "gamma", "lognormal", "normal", or "uniform".

• p1: The first hyperparameter:

  • "gamma": shape

  • "lognormal": meanlog (mean on the log scale)

  • "normal": mean

  • "uniform": lower bound

• p2: The second hyperparameter:

  • "gamma": rate

  • "lognormal": sdlog (standard deviation on the log scale)

  • "normal": sd (standard deviation)

  • "uniform": upper bound

target

Character vector ("gdata" or "ldata") that determines if the ABC-SMC method estimates parameters in model@gdata (global data) or in model@ldata (local data).

pars

An integer vector with the indices of the parameters in target that are being estimated.

nprop

An integer vector with the number of simulated proposals (particles) generated in each generation.

fn

A function used to calculate summary statistics from the simulated trajectory and compute the distance for each particle. See abc for details on the required function signature.

tolerance

A numeric matrix (number of summary statistics $\times$ number of generations) where each column contains the tolerances for a generation and each row contains a sequence of gradually decreasing tolerances.

x

A numeric array (number of particles $\times$ number of parameters $\times$ number of generations) with the parameter values for the accepted particles in each generation. Each row is one particle.

weight

A numeric matrix (number of particles $\times$ number of generations) with the weights for the particles x in the corresponding generation.

distance

A numeric array (number of particles $\times$ number of summary statistics $\times$ number of generations) with the distance for the particles x in each generation. Each row contains the distance for a particle and each column contains the distance for a summary statistic.

ess

A numeric vector with the effective sample size (ESS) in each generation. The effective sample size is computed as

$\left(\sum_{i=1}^N\!(w_{g}^{(i)})^2\right)^{-1},$

where $w_{g}^{(i)}$ is the normalized weight of particle $i$ in generation $g$.

init_model

An optional function that, if non-NULL, is applied before running each proposal. The function must accept one argument of type SimInf_model with the current model of the fitting process and return a modified model. This function can be useful to specify the initial state of u0 or v0 of the model before running a trajectory with proposed parameters.

See Also

abc for the main ABC function and continue_abc for continuing an ABC run.

Create a SimInf_events object

Description

The argument events must be a data.frame with the following columns:

event

Four event types are supported by the current solvers: exit, enter, internal transfer, and external transfer. When assigning the events, they can either be coded as a numerical value or a character string: exit; 0 or 'exit', enter; 1 or 'enter', internal transfer; 2 or 'intTrans', and external transfer; 3 or 'extTrans'. Internally in SimInf, the event type is coded as a numerical value.

time

When the event occurs i.e., the event is processed when time is reached in the simulation. Can be either an integer or a Date vector. A Date vector is coerced to a numeric vector as days, where t0 determines the offset to match the time of the events to the model tspan vector.

node

The node that the event operates on. Also the source node for an external transfer event. 1 <= node[i] <= Number of nodes.

dest

The destination node for an external transfer event i.e., individuals are moved from node to dest, where 1 <= dest[i] <= Number of nodes. Set event = 0 for the other event types. dest is an integer vector.

n

The number of individuals affected by the event. n[i] >= 0.

proportion

If n[i] equals zero, the number of individuals affected by event[i] is calculated by sampling the number of individuals from a binomial distribution using the proportion[i] and the number of individuals in the compartments. Numeric vector. 0 <= proportion[i] <= 1.

select

To process an event[i], the compartments affected by the event are specified with select[i] together with the matrix E, where select[i] determines which column in E to use. The specific individuals affected by the event are sampled from the compartments corresponding to the non-zero entries in the specified column in E[, select[i]], where select is an integer vector.

shift

Determines how individuals in enter, internal transfer and external transfer events are shifted to enter another compartment. The sampled individuals are shifted according to column shift[i] in matrix N i.e., N[, shift[i]], where shift is an integer vector. See above for a description of N. Unsued for the other event types.

Usage

SimInf_events(E = NULL, N = NULL, events = NULL, t0 = NULL)

Arguments

Value

S4 class SimInf_events

Examples

## Let us illustrate how movement events can be used to transfer
## individuals from one node to another.  Use the built-in SIR
## model and start with 2 nodes where all individuals are in the
## first node (100 per compartment).
u0 <- data.frame(S = c(100, 0), I = c(100, 0), R = c(100, 0))

## Then create 300 movement events to transfer all individuals,
## one per day, from the first node to the second node. Use the
## fourth column in the select matrix where all compartments
## can be sampled with equal weight.
events <- data.frame(event      = rep("extTrans", 300),
                     time       = 1:300,
                     node       = 1,
                     dest       = 2,
                     n          = 1,
                     proportion = 0,
                     select     = 4,
                     shift      = 0)

## Create an SIR model without disease transmission to
## demonstrate the events.
model <- SIR(u0      = u0,
             tspan  = 1:300,
             events = events,
             beta   = 0,
             gamma  = 0)

## Run the model and plot the number of individuals in
## the second node.  As can be seen in the figure, all
## indivuduals have been moved to the second node when
## t = 300.
plot(run(model), index = 1:2, range = FALSE)

## Let us now double the weight to sample from the 'I'
## compartment and rerun the model.
model@events@E[2, 4] <- 2
plot(run(model), index = 1:2, range = FALSE)

## And much larger weight to sample from the I compartment.
model@events@E[2, 4] <- 10
plot(run(model), index = 1:2, range = FALSE)

## Increase the weight for the R compartment.
model@events@E[3, 4] <- 4
plot(run(model), index = 1:2, range = FALSE)

Description

Class to hold data for scheduled events that modify the discrete state of individuals in a node at a pre-defined time t.