This vignette covers how to work with network objects, edgelists,
multi-layer models, clique and observed network layers, and partnership
histories in EpiModel network models with custom extension modules. It
assumes familiarity with setting up and running network models with
netsim() and with the extension API. See the Network Modeling for
Epidemics (NME) course materials and the EpiModel Gallery
for background.
For working with nodal attributes and epidemic summary statistics, see the companion vignette Working with Custom Attributes and Summary Statistics in EpiModel.
EpiModel supports two storage modes for networks, controlled by the
tergmLite parameter in control.net():
Full mode (tergmLite = FALSE, the
default): Networks are stored as networkDynamic objects,
which preserve the complete history of edge activations and
deactivations. This allows extraction of the full dynamic network after
simulation. However, networkDynamic objects consume
substantial memory.
tergmLite mode (tergmLite = TRUE):
Networks are stored as lightweight networkLite objects
containing only the current edgelist and nodal attributes. This provides
a 20–50x performance improvement and much lower memory usage, making it
essential for large-scale research models. The trade-off is that
get_network() returns a networkLite (a
snapshot) rather than a full dynamic network history.
Most extension modules work identically under both modes because they access networks through EpiModel’s accessor functions rather than manipulating network objects directly.
Inside a custom module, use the get_network() and
set_network() accessors to work with network objects:
# Get the network for layer 1
nw <- get_network(dat, network = 1)
# After modifying a network, set it back
dat <- set_network(dat, network = 1, nw = nw)These accessors handle the tergmLite vs. full-mode distinction
internally, so your module code works under both storage modes. In full
mode, get_network() returns a networkDynamic
object; in tergmLite mode, it returns a networkLite.
In practice, you rarely need to access network objects directly. Instead, use the edgelist accessor functions described below, which also work correctly under both storage modes.
After a netsim call, extract network objects with
get_network():
sim <- netsim(est, param, init, control)
# Extract the network from simulation 1, network layer 1
nw <- get_network(sim, sim = 1, network = 1)
# Collapse to a static cross-section at time step 50 (full mode only)
nw_at_50 <- get_network(sim, sim = 1, collapse = TRUE, at = 50)In full mode, get_network() returns a
networkDynamic object. In tergmLite mode, it returns a
networkLite object representing the final state. The
collapse and at arguments are only available
in full mode.
Note: Network objects are only saved in the output
when save.network is TRUE in
control.net() (the default for full mode). In tergmLite
mode, there is no network history to save.
The transmission matrix records every transmission event during the simulation:
transmat <- get_transmat(sim, sim = 1)This returns a data.frame with columns including
at (time step), sus (ID of the newly infected
node), inf (ID of the infecting node), network
(the layer on which the transmission occurred), infDur
(duration of infector’s infection), transProb,
actRate, and finalProb. Transmission matrices
are saved by default (save.transmat = TRUE in
control.net()).
Current edgelists are the set of active partnerships at the present time step. These are the most commonly used network data structures inside extension modules.
get_edgelist() returns the current edgelist for a given
network as a two-column matrix of positional IDs:
el <- get_edgelist(dat, network = 1)Each row is an active partnership. Column 1 is the positional ID of the “head” node; column 2 is the “tail” node. This function works identically under both storage modes.
get_edgelists_df() combines edgelists from multiple
network layers into a single data.frame with a
network column:
# All networks
el_all <- get_edgelists_df(dat, networks = NULL)
# Specific networks
el_12 <- get_edgelists_df(dat, networks = c(1, 2))The returned data.frame has columns head,
tail, and network.
The discordant edgelist identifies partnerships where partners have different values of a status attribute—the key data structure for modeling transmission. For example, in an SI model, discordant edges are those where one partner is susceptible and the other is infected:
disc_el <- get_discordant_edgelist(
dat,
status.attr = "status",
head.status = "i",
tail.status = "s"
)The returned data.frame has columns head,
tail, head_status, tail_status,
and network. Both orderings are captured: if node A
(infected) is partnered with node B (susceptible), the edge appears
regardless of which is the “head” vs “tail” in the underlying
network.
See also discord_edgelist() for the original, simpler
version of this function used in built-in models.
EpiModel uses two ways to reference nodes:
By position: Think of it like a row number in a
spreadsheet. get_attr(dat, "active", posit_ids = 3)
accesses the third node’s value directly. This is the standard way to
look up node information and is very fast. In a model with 100 nodes,
positions range from 1 to 100. When nodes depart, they may be dropped
from the vectors, freeing their position for new arrivals.
By unique_id: A globally unique
integer attribute assigned to each node at creation and never reused.
Slower to look up, but allows referencing nodes that have already
departed. Used by cumulative edgelists and attribute histories.
Conversion between the two systems is handled internally by EpiModel.
The get_unique_ids() and get_posit_ids()
functions perform the conversion. See
help("unique_id-tools", package = "EpiModel") for
details.
A multi-layer model has several edge sets over one node set, such as
main and casual sexual partnerships, sexual and needle-sharing
partnerships, or household and community contacts. Each layer estimated
with netest() has its own formation and dissolution model,
and infection can be transmitted across the edges of any layer. The node
set, including the network size and the nodal attributes, is shared by
all layers. A layer may also be a clique layer built with
netclique() or an observed network wrapped with
netcensus(); both are described in the sections that
follow. The Multi-Layer
Networks chapter of the NME course materials works through a
complete two-layer model; this section summarizes the mechanics.
Each layer is estimated with its own netest() call on
the same starting network, so that every layer carries the same nodes
and nodal attributes. The fitted layers are then passed to
netsim() as a list:
nw <- network_initialize(n = 1000)
nw <- set_vertex_attribute(nw, "race", rep(0:1, each = 500))
est_main <- netest(nw, formation = ~edges + nodematch("race"),
target.stats = c(300, 240),
coef.diss = dissolution_coefs(~offset(edges), duration = 100))
est_casl <- netest(nw, formation = ~edges + nodematch("race"),
target.stats = c(200, 150),
coef.diss = dissolution_coefs(~offset(edges), duration = 10))
sim <- netsim(list(est_main, est_casl), param, init, control)The order of the list sets the network index of each layer: here the
main layer is network 1 and the casual layer is network 2. Every
accessor below refers to a layer by this index. netsim()
reads the nodal attributes from the network of the first layer in the
list. Each layer is diagnosed separately with netdx(), as
in a single-layer model.
The multilayer() function specifies one value per layer,
in the order of the layer list. Four control.net()
arguments accept it: nwstats.formula (the network
statistics recorded for each layer), set.control.ergm,
set.control.tergm, and
tergmLite.track.duration. In param.net(),
inf.prob, inf.prob.g2, and
act.rate accept it, and the built-in infection modules then
apply each layer’s value to the discordant edges of that layer. An
argument that is not a multilayer object applies to every
layer:
param <- param.net(inf.prob = multilayer(0.2, 0.1),
act.rate = multilayer(2, 1))
control <- control.net(type = "SI", nsteps = 100, nsims = 1,
tergmLite = TRUE, resimulate.network = TRUE,
nwstats.formula = multilayer("formation",
~edges + degree(0:3)))Each entry of a multilayer parameter may itself vary
with the duration of infection, given as a
by_infection_duration() object:
inf.prob = multilayer(by_infection_duration(c(0.5, 0.2)), 0.1)
gives the first layer a probability of 0.5 in the first time step of an
infection and 0.2 after.
In a custom infection module, get_param() returns the
multilayer object, a list with one element per layer, and
the network column of the discordant edgelist selects the
element that applies to each edge:
inf.prob <- get_param(dat, "inf.prob")
del <- get_discordant_edgelist(dat, status.attr = "status",
head.status = "i", tail.status = "s")
del$transProb <- unlist(inf.prob)[del$network] # one fixed value per layerThe network accessors take the layer index in their
network argument:
get_network(dat, network = 2) and
get_edgelist(dat, network = 2) read one layer, while
get_edgelists_df() and
get_discordant_edgelist() combine the layers and add a
network column (see Current Edgelists above).
dat$num.nw holds the number of layers.get_network(sim, network = 2) extracts one layer, and
get_nwstats(sim, network = 2),
print(sim, network = 2), and
plot(sim, type = "formation", network = 2) report its
network statistics.network column, so
table(get_transmat(sim)$network) counts the transmissions
on each layer.Layers estimated separately are independent: the edges of one layer do not affect the formation of edges in another. Dependence between layers enters through a nodal attribute that summarizes one layer and appears in the formation model of another. Here, the formation of casual partnerships depends on whether a node has a main partner:
# nw carries a starting value of deg.main (1 if the node has a main partner)
est_casl <- netest(nw, formation = ~edges + nodematch("race") +
nodefactor("deg.main"),
target.stats = c(200, 150, 40),
coef.diss = dissolution_coefs(~offset(edges), duration = 10))
update_deg_main <- function(dat, at, network) {
if (network == 1) {
deg <- get_degree(get_edgelist(dat, network = 1))
dat <- set_attr(dat, "deg.main", pmin(deg, 1))
}
return(dat)
}
control <- control.net(type = "SI", nsteps = 100, nsims = 1,
tergmLite = TRUE, resimulate.network = TRUE,
dat.updates = update_deg_main)The starting values of the attribute should be consistent with the
main layer, and they must be on the network of the first layer in the
list, from which netsim() reads nodal attributes. The
attribute must then be kept current as the layers are resimulated. The
dat.updates argument of control.net() takes a
function of dat, at, and network
that netsim() calls at every time step, including
initialization, before the first layer is resimulated
(network = 0) and after each layer is resimulated
(network = 1, 2, ...). The function above recomputes
deg.main after the main layer is redrawn, so that the
casual layer is resimulated against the current main partnerships. When
each layer depends on the other, the function also recomputes the
attribute that summarizes the second layer; the NME chapter covers that
case, including how to seed the starting values with
san().
This pattern assumes tergmLite = TRUE, under which the
resimulation reads nodal attributes from dat. In full mode,
each network object carries its own copy of the nodal attributes, which
nwupdate.net() refreshes from dat once per
time step, so an attribute set in dat.updates does not
reach the next layer’s resimulation in the same time step. Because main
partnerships begin and end while casual partnerships formed under the
earlier value of deg.main persist, the cross-layer
statistic can drift from its target during the simulation; recording it
through nwstats.formula shows by how much.
Some contact structures are groups by definition: every pair of
people in a household, classroom, hospital ward, or ship cabin is in
contact for as long as they share it. netclique() builds a
network layer of this kind from a grouping attribute, connecting every
pair of nodes that share a value so that each group is a clique.
netsim() accepts the result anywhere in its list of layers,
alongside layers estimated with netest(), and the edges of
the layer never form or dissolve on their own.
A clique layer is an addition to netest(), not a
substitute for it. An ERGM remains the model for any network whose ties
depend on nodal and dyadic predictors (degree, mixing by attribute,
clustering), whether those ties turn over quickly, slowly, or not at
all; a long partnership duration in dissolution_coefs()
keeps a netest layer close to fixed. A clique layer has no
tie-formation process to estimate. An ERGM could reproduce the cliques
only through a nodematch term on the grouping attribute
targeted at its maximum, where no finite coefficient exists. See
help("netclique", package = "EpiModel") for the full
rationale.
The grouping attribute is a vertex attribute (integer, numeric, or
character, but not a factor) whose values partition the nodes into
groups. Nodes with a missing (NA) value belong to no group
and are isolates on the layer. Two helpers build the attribute:
sample_groups() builds a population one group at a time
from a table of group types, such as household compositions by age
group, and returns each node’s group ID together with the attributes its
group type implies.assign_groups() solves the reverse problem. It assigns
group IDs to a population whose attributes already exist, drawing group
sizes from a target distribution and keeping dependent members (such as
children) in groups with an anchor member (such as an adult).Here, a population of 1,000 is drawn from four household types, and the same starting network is used for a household clique layer and a community TERGM layer:
hh <- sample_groups(1000, c("adult" = 0.25, "adult adult" = 0.35,
"adult adult child" = 0.25,
"adult adult child child" = 0.15),
attr.name = "age")
nw <- network_initialize(n = 1000)
nw <- set_vertex_attribute(nw, "age", hh$age)
nw <- set_vertex_attribute(nw, "hh_id", hh$group)
est_hh <- netclique(nw, group.attr = "hh_id", arrivals = "join")
est_com <- netest(nw, formation = ~edges + nodematch("age"),
target.stats = c(400, 300),
coef.diss = dissolution_coefs(~offset(edges), duration = 20))
param <- param.net(inf.prob = multilayer(0.3, 0.05), act.rate = 1)
sim <- netsim(list(est_hh, est_com), param, init, control)The household layer is network 1 and the community layer is network
2, and the multilayer() parameter gives household contacts
the higher per-act transmission probability. Since netsim()
reads nodal attributes from the first layer in the list, the grouping
attribute should also be set on the network passed to
netest() for the other layers, as it is here; when it is
not, netsim() copies it from the clique layer.
Printing a netclique object shows the number of groups,
the edge count, the mean degree, and the group size distribution, and
print(est_hh, by = "age") adds the mean degree by a nodal
attribute. There is no model to diagnose, so netdx() does
not accept a netclique object.
A clique layer is stored like any other layer: as an edgelist in
tergmLite mode and as a networkDynamic object in full mode.
The accessors in this vignette (get_network(),
get_edgelist(), get_edgelists_df(),
get_discordant_edgelist()) read it without special
handling, and the built-in infection modules transmit across its edges
as on any other layer. The network resimulation and the edges correction
skip the layer, so its edges change only when nodes depart, when nodes
arrive, and when a module calls move_to_group(). With the
default nwstats.formula = "formation", the network
statistic recorded for the layer is its edge count.
Departing nodes are removed from a clique layer together with their
edges, as on every other layer. Arriving nodes are placed on the layer
under the rule chosen with the arrivals argument of
netclique(). The rule is applied within
arrive_nodes(), called from the built-in
nwupdate.net() module after the arrivals module has created
the new nodes and set their other attributes:
"isolate" (the default): the new node has no edges on
the layer and its grouping attribute is NA. A custom module
may place it later with move_to_group()."new": each new node starts a group of size one, with a
fresh group ID."join": each new node joins an existing group and is
connected to every active member of that group, including other nodes
joining the same group in the same time step. By default, the group is
drawn with probability proportional to its current size.Under "join", an arrivals.FUN function
replaces the default draw. It is called as
arrivals.FUN(dat, at, new_ids, network), where
new_ids holds the positional IDs of the new nodes, and it
returns one group ID per new node (NA leaves that node as
an isolate). Because it receives dat, it can read any nodal
attribute. Here, each newborn joins a household that already has a
child:
est_hh <- netclique(nw, group.attr = "hh_id", arrivals = "join",
arrivals.FUN = function(dat, at, new_ids, network) {
active <- get_attr(dat, "active")
hh <- get_attr(dat, "hh_id")
age <- get_attr(dat, "age")
existing <- active == 1
existing[new_ids] <- FALSE
pool <- hh[which(existing & age == "child" & !is.na(hh))]
pool[sample.int(length(pool), length(new_ids), replace = TRUE)]
})The pool leaves out new_ids. When the arrivals module
gives new nodes a default value of the grouping attribute, such as
0, a pool taken from all nodes would contain that value,
and every arrival that drew it would join the other arrivals in a
spurious group. Define arrivals.FUN at the top level of a
script or in a package: the layer keeps the function with its enclosing
environment, and so does every simulation run from it.
A custom arrivals module does not need to set the grouping attribute:
the new nodes receive NA unless the module gives the
attribute a default value, and the layer’s rule then places them either
way. A module that does set the attribute itself should be paired with
arrivals = "join" and an arrivals.FUN that
returns those values,
function(dat, at, new_ids, network) get_attr(dat, "hh_id")[new_ids],
so that the new nodes are connected to the members of their groups.
The clique edges are built from the grouping attribute once, when the
simulation starts. Setting the attribute with set_attr()
afterward does not rewire the layer: the old edges remain and the layer
no longer matches the groups. move_to_group() sets the
attribute and the edges together, in either storage mode. Each moved
node loses its edges to the members of its old group and is connected to
every active member of its new group. In a model where age
is in years, a module in which young adults leave home to start
households of their own is:
leave_home <- function(dat, at) {
active <- get_attr(dat, "active")
age <- get_attr(dat, "age")
hh_id <- get_attr(dat, "hh_id")
elig <- which(active == 1 & age >= 18 & age < 30)
movers <- elig[runif(length(elig)) < 0.01]
if (length(movers) > 0) {
new_ids <- max(hh_id, na.rm = TRUE) + seq_along(movers)
dat <- move_to_group(dat, ids = movers, group = new_ids)
}
return(dat)
}The ids argument takes positional IDs. The
group argument takes one group ID per node, or a single
value for all of them: an ID in use joins that group, an unused ID
starts a new group, and NA takes the node out of its group.
In a model with more than one clique layer, the network
argument names the layer to change.
The cumulative edgelist (see below) records the initial clique edges
with a start of 0. Edges added by arrivals and moves, and
edges ended by moves and departures, are recorded when the network
resimulation module next updates the cumulative edgelist, which it does
for every layer at every time step. A module that reads the cumulative
edgelist later in the same time step as a move should call
update_cumulative_edgelist() after the move.
netest() turns partial, usually egocentric, network data
into a generative model that netsim() can simulate from.
When the whole network has been observed (every node and every contact,
and for a dynamic network every time step), there is nothing to
estimate: the observed object is what netsim() would
otherwise have to generate. netcensus() wraps such a
network, either a static network object or a
networkDynamic object with edge spells, as a layer that
netsim() accepts anywhere in its list of layers. Sensor,
proximity-logger, contact-tracing, and animal-tracking datasets are the
usual sources.
Like a clique layer, an observed layer is an addition to
netest(), not a substitute for it. It is the wrong tool for
a sample from which one wants to generalize to a population, and it is
not a way to model how the observed ties arise; a network whose ties are
to be reproduced from their predictors is estimated with
netest(), whatever their turnover. The Epidemics
over Observed Networks chapter of the NME course materials works
through a complete example.
library(networkDynamicData)
data(concurrencyComparisonNets)
obs <- netcensus(base)
obsnetcensus() carries the vertex attributes of the
observed network into the simulation as nodal attributes, wherever the
layer sits in the list of layers; an attribute that the first layer also
holds must have the same values on both. The network must be undirected
and not bipartite, without loops or parallel edges (a pair of nodes in
contact over several periods has one edge with several activity spells).
Temporally extended vertex attributes (those ending in
.active), such as the status.active attribute
stored on base, are dropped with a message, since
netsim() keeps its own disease status. Printing a dynamic
census shows the observation window and the number of edges active per
step across it; printing a static census shows its edge count, mean
degree, and number of isolates. As with a clique layer, there is no
model to diagnose, so netdx() does not accept a
netcensus object.
A static census has the same edges at every time step. For a dynamic
census, EpiModel time step at reads the edges active at
time at of the observed object, so the simulation clock is
the observation clock. A census recorded in other units, such as seconds
or calendar dates, needs its spell times rescaled to time steps first.
The window argument of netcensus() sets the
observation window, read like a networkDynamic spell as the
half-open interval [start, end): c(1, 11)
covers time steps 1 to 10. It defaults to the
net.obs.period network attribute when the object has one
with finite bounds, and otherwise to the range of the finite edge spell
times; for base, the window is [2, 102).
Past the end of the window, the layer holds only the edges whose
spells extend beyond it. Spells censored at the end of the observation
(a terminus of Inf, as in tergm and EpiModel
output such as base) stay active, so the edge set is
frozen; spells that end with the observation, as in most contact data,
leave the layer empty. Keep the time steps of control.net()
within the window; netsim() warns when the simulation ends
before the window begins or nsteps reaches its end:
param <- param.net(inf.prob = 0.5, act.rate = 1)
init <- init.net(i.num = 10)
control <- control.net(type = "SI", nsteps = 100, nsims = 5,
resimulate.network = FALSE)
sim <- netsim(obs, param, init, control)An observed layer is skipped by the network resimulation whatever the
value of resimulate.network. Setting it to
FALSE is enough here because the only layer is observed;
the setting matters only when an observed layer sits next to an
estimated one.
Without tergmLite, the observed
networkDynamic object is used as stored, and
get_edgelist() reads the edges active at the current time
step from it. Under tergmLite, the edgelist of the layer is
replaced at every time step with the edges active in the observed
object. Either way, the accessors read the layer as any other, and the
built-in infection modules transmit across its edges. Next to other
layers, the per-layer inf.prob and act.rate
are set through multilayer():
# nw_obs: an observed static network on the same nodes as the estimated layer
sim <- netsim(list(netcensus(nw_obs), est),
param.net(inf.prob = multilayer(0.3, 0.1), act.rate = 1),
init, control)An observed layer differs from an estimated one in a few other respects:
netsim() refuses vital dynamics in a model with an observed
layer, and arrive_nodes() and depart_nodes()
stop if a custom module tries to add or remove nodes. Vertex activity
spells are not used: every node is present throughout, and a node absent
from part of the observation simply has no contacts then.nwstats.formula = "formation", the statistic recorded for
the layer is its edge count, summarized at every time step for a dynamic
census. There are no target statistics, so print(sim) and
plot(sim, type = "formation") show the recorded values
without targets.tergmLite (tergmLite.track.duration = TRUE) is
refused for a dynamic census, whose observed spells already carry the
edge durations.get_nwparam(dat, network = 1)$census.nw.tergmLite, but the current step’s edges
in full mode, which reads the networkDynamic object
directly. Modules added to control.net() under new names
run before the built-in modules unless module.order places
them later.The cumulative edgelist is a historical record of all edges in a network, including the time steps when each edge started and stopped. This allows querying both current and past partnerships—essential for contact tracing, partnership duration analysis, and reachability analysis.
The cumulative edgelist follows a four-step lifecycle. The same data structure is produced once and read at different stages:
control.net(): set
cumulative.edgelist = TRUE (in-memory tracking during the
run) and, if you want to use the result after netsim()
returns, also set save.cumulative.edgelist = TRUE.control$truncate.el.cuml
to decide how much dissolved-edge history to retain. Custom modules that
mutate the network outside the TERGM machinery should also call
update_cumulative_edgelist() after the change.get_cumulative_edgelist() /
get_cumulative_edgelists_df(), or with one of the derived
helpers below (get_partners(),
get_cumulative_degree()).sim$cumulative.edgelist[[s]] directly on the returned
netsim object.The control settings, accessors, and helpers all share the same
column convention: head/tail (or
index/partner) by unique ID,
start/stop time steps (inclusive, with
NA stop for active edges), and
network index for multi-layer models.
Cumulative edgelist tracking must be explicitly enabled in
control.net():
control <- control.net(
type = "SI",
nsims = 1,
nsteps = 100,
cumulative.edgelist = TRUE, # Enable in-memory tracking during the run
truncate.el.cuml = 0, # Drop dissolved edges immediately (default)
save.cumulative.edgelist = TRUE, # Attach the result to the returned netsim
verbose = FALSE
)Without cumulative.edgelist = TRUE, calls to
update_cumulative_edgelist() silently do nothing, and calls
to get_cumulative_edgelist() raise an error.
The truncate.el.cuml control sets the default truncation
passed to every automatic update_cumulative_edgelist()
call. The default (0) keeps only currently active edges,
which is enough for tracking active-edge start times while keeping
memory low. Use Inf to keep full history, or a positive
integer to retain dissolved edges for that many steps after they
ended.
The cumulative edgelist must be updated at each time step after network resimulation. In a custom module or at the end of initialization:
dat <- update_cumulative_edgelist(dat, network = 1, truncate = Inf)The truncate argument controls memory usage:
truncate = Inf: Keep the full history of all edges (no
removal).truncate = 0 (the default): Keep only currently active
edges. Use this if you only need to track active edge start times.truncate = N: Remove edges that ended more than
N time steps ago. This balances historical depth with
memory use.To update all networks in a multi-layer model:
for (n_network in seq_len(dat$num.nw)) {
dat <- update_cumulative_edgelist(dat, n_network, truncate = 100)
}Inside a custom module, where dat is the live
netsim_dat object:
el_cuml <- get_cumulative_edgelist(dat, network = 1)The returned tibble has four columns:
head: the unique_id of the first
node.tail: the unique_id of the second
node.start: the time step when the edge formed.stop: the last time step the edge was active, or
NA if the edge is still active.An edge with start = 5 and stop = 12
existed from steps 5 through 12, inclusive. When no edges have been
recorded yet, the function still returns a tibble with
these four columns and zero rows. If
cumulative.edgelist = FALSE, the function raises an error
rather than silently returning empty.
el_cumls <- get_cumulative_edgelists_df(dat, networks = NULL)The networks argument accepts a vector of network
indices or NULL (all networks). The returned
data.frame adds a network column identifying
which network layer each edge belongs to.
When save.cumulative.edgelist = TRUE, the cumulative
edgelist is attached to the netsim return object as
sim$cumulative.edgelist, a list with one element per
simulation:
el_cuml <- sim$cumulative.edgelist[[1]]Each element is the data.frame produced by
get_cumulative_edgelists_df() (head,
tail, start, stop,
network). Do not call the accessor functions on a processed
netsim object: the live dat$run state they
depend on is dropped during process_out.net(). Note also
that truncate_sim() keeps the edges of
sim$cumulative.edgelist that ended before at.
With reset.time = TRUE (the default), their
start and stop are shifted to the new time
steps along with the rest of the output, so
is.na(stop) | stop >= 1 keeps the edges active in the
kept window; with reset.time = FALSE the times are
unchanged and the filter is
is.na(stop) | stop >= at.
get_partners() extracts the partners of specified nodes
from the cumulative edgelist:
partner_list <- get_partners(
dat,
index_posit_ids,
networks = NULL,
truncate = Inf,
only.active.nodes = FALSE
)Arguments:
dat: the main list object.index_posit_ids: a vector of positional IDs for the
nodes of interest (the “indexes”).networks: which network layers to search
(NULL for all).truncate: only include edges that ended within this
many steps (filters by edge age).only.active.nodes: if TRUE, exclude
partnerships with inactive (departed) nodes.The output is similar to get_cumulative_edgelists_df()
but with columns index and partner (both
containing unique IDs) instead of head and
tail. Note that indexes are specified by positional ID but
the output uses unique IDs, since partners may include nodes that have
already departed.
get_cumulative_degree() counts the number of distinct
partners each node has had over the tracked history:
cum_degree <- get_cumulative_degree(
dat,
index_posit_ids = 1:50,
networks = NULL,
truncate = Inf,
only.active.nodes = FALSE
)This returns a data.frame with columns
index_pid (positional ID) and degree
(cumulative partner count). It wraps get_partners() and
counts unique partners per index.
Reachability functions determine which nodes can be connected through chains of partnerships over a time window. These are useful for outbreak investigation (forward reachability: who could a node have infected?) and source tracing (backward reachability: who could have infected a node?).
These functions operate on cumulative edgelist objects directly—not
on dat—and are typically used for post-simulation analysis.
Most often, the input comes from the saved
sim$cumulative.edgelist[[s]] slot of a netsim
run; see the Converting External networkDynamic
Objects and Deduplicating Across Sources subsections below
for two ways to construct input from other starting points.
el_cuml <- get_cumulative_edgelist(dat, network = 1)
fwd <- get_forward_reachable(
el_cuml,
from_step = 1,
to_step = 52,
nodes = c(10, 25, 42), # NULL for all nodes with edges
dense_optim = "auto"
)Returns a list with two elements:
reached: a named list where each element contains the
set of nodes reachable from each index node through chains of
partnerships active during [from_step, to_step].lengths: a matrix with one row per node and one column
per time step (plus an initial column), showing how the reachable set
grows over time. This allows back-calculating the distance to each
reachable node.Nodes are identified by unique ID in the output (named as
node_ID). Nodes with no edges during the analysis period
are excluded from the output; their forward reachable set is just
themselves (size 1).
bkwd <- get_backward_reachable(
el_cuml,
from_step = 1,
to_step = 52,
nodes = c(10, 25, 42)
)Same interface and output structure as
get_forward_reachable(), but follows partnerships backward
in time. This answers the question: which nodes could have reached this
node through a chain of partnerships during the specified period?
networkDynamic ObjectsIf you have a networkDynamic object that was produced
outside of a cumulative.edgelist = TRUE run (for example,
from a non-tergmLite simulation, from a saved RDS file, or
built manually), as_cumulative_edgelist() converts it into
the same data.frame shape the reachability functions
expect:
nd <- get_network(sim, sim = 1, network = 1) # full-mode netsim only
el_cuml <- as_cumulative_edgelist(nd)
fwd <- get_forward_reachable(el_cuml, from_step = 1, to_step = 100)This is the recommended entry point for post-hoc reachability
analysis on simulations that did not enable cumulative-edgelist tracking
up front, provided the run kept the full networkDynamic
history (i.e., not tergmLite).
The reachability functions assume non-overlapping spells per
(head, tail) pair. When concatenating cumulative edgelists
from multiple sources (e.g., separate simulation segments, or several
as_cumulative_edgelist() conversions stitched together),
use dedup_cumulative_edgelist() to merge overlapping spells
into a single row:
el_cuml <- dedup_cumulative_edgelist(rbind(el_run1, el_run2))The reachability functions are efficient for large networks because
they operate on cumulative edgelists (much smaller than full
networkDynamic objects). For sets of more than 5 nodes,
they are faster than iterating over tsna::tPath(). The
dense_optim argument controls an adjacency-list
optimization that helps with dense networks; "auto" enables
it when the number of edges exceeds the number of nodes.