## ----echo=FALSE---------------------------------------------------------------
vers <- packageVersion("EpiModel")
year <- format(Sys.time(), "%Y")

## ----step1, results = "hide", message = FALSE---------------------------------
library(EpiModel)
set.seed(12345)

nw <- network_initialize(n = 500)

formation <- ~edges + concurrent + degrange(from = 4)
target.stats <- c(175, 110, 0)
coef.diss <- dissolution_coefs(
  dissolution = ~offset(edges),
  duration = 50   # mean partnership duration in time steps
)

est <- netest(nw, formation, target.stats, coef.diss, verbose = FALSE)

## ----step2, results = "hide"--------------------------------------------------
dx <- netdx(est, nsims = 5, nsteps = 500, verbose = FALSE)

## ----step2_plot1, fig.width = 8-----------------------------------------------
plot(dx)

## ----step2_plot2--------------------------------------------------------------
plot(dx, type = "duration")

## ----step3, results = "hide"--------------------------------------------------
param <- param.net(inf.prob = 0.4, act.rate = 2, rec.rate = 0.05)
init <- init.net(i.num = 10)
control <- control.net(type = "SIS", nsims = 5, nsteps = 500, verbose = FALSE)

sim <- netsim(est, param, init, control)

## ----step4_epi, fig.width = 8-------------------------------------------------
plot(sim, sim.lines = TRUE, sim.alpha = 0.3,
     mean.lwd = 3, qnts = 0.5, main = "SIS Epidemic on a Dynamic Network")

## ----step4_prev---------------------------------------------------------------
plot(sim, y = "i.num", popfrac = TRUE, sim.lines = TRUE, sim.alpha = 0.3,
     mean.lwd = 3, qnts = 0.5, main = "Infection Prevalence",
     ylab = "Prevalence", legend = FALSE)

## ----step4_summary------------------------------------------------------------
summary(sim, at = 500)

## ----step4_df-----------------------------------------------------------------
d <- as.data.frame(sim)
head(d)

## ----step4_transmat-----------------------------------------------------------
# Run a separate SI model seeded with a single infection for a cleaner tree
set.seed(42)
param2 <- param.net(inf.prob = 0.5, act.rate = 2)
init2 <- init.net(i.num = 5)
control2 <- control.net(type = "SI", nsims = 1, nsteps = 60, verbose = FALSE)
sim2 <- netsim(est, param2, init2, control2)

tm <- get_transmat(sim2)
head(tm)

## ----step4_tree, fig.width = 8, fig.height = 5--------------------------------
tmPhylo <- as.phylo.transmat(tm)
par(mar = c(2, 1, 2, 1))
plot(tmPhylo, show.node.label = TRUE, root.edge = TRUE, cex = 0.5,
     main = "Transmission Tree")

## ----extension, eval = FALSE--------------------------------------------------
# # Example: a custom progression module for an SEIR model
# progress_module <- function(dat, at) {
#   active <- get_attr(dat, "active")
#   status <- get_attr(dat, "status")
# 
#   # E -> I transition
#   ids.EtoI <- which(active == 1 & status == "e" & runif(length(status)) < 0.1)
#   status[ids.EtoI] <- "i"
# 
#   # I -> R transition
#   ids.ItoR <- which(active == 1 & status == "i" & runif(length(status)) < 0.02)
#   status[ids.ItoR] <- "r"
# 
#   dat <- set_attr(dat, "status", status)
#   dat <- set_epi(dat, "ei.flow", at, length(ids.EtoI))
#   dat <- set_epi(dat, "ir.flow", at, length(ids.ItoR))
#   return(dat)
# }

## ----eval=FALSE---------------------------------------------------------------
# help(package = "EpiModel")