make_data: Simulate time series data

Description Usage Arguments Value Model Randomness Examples

View source: R/make_data.R

Description

make_data() simulates time series data according to a system of SIR equations and a supplied list of parameter values. Observations are recorded at equally spaced time points tk = t0+kΔt (for k = 0,...,n), where Δt denotes the observation interval.

Usage

1
2
3
4
5
6
7
make_data(
  par_list = list(),
  n = 20 * 365/7,
  with_dem_stoch = TRUE,
  seed = NA,
  ode_control = list(method = "lsoda", rtol = 1e-06, atol = 1e-06)
)

Arguments

par_list

A list containing:

dt_weeks

[ Δt ] Observation interval in weeks.

t0

[ t0 ] Time of the first observation in units Δt.

prep

[ prep ] Case reporting probability.

trep

[ trep ] Case reporting delay in units Δt.

hatN0

[ Ñ0 ] Population size at time t = 0 years.

N0

[ N0 ] Population size at time t = t0.

S0

[ S0 ] Number of susceptibles at time t = t0.

I0

[ I0 ] Number of infecteds at time t = t0.

nu

[ νc ] Birth rate expressed per unit Δt and relative to Ñ0.

mu

[ μc ] Natural mortality rate expressed per unit Δt and per capita.

tgen

[ tgen ] Mean generation interval of the disease of interest in units Δt.

beta_mean

[ ⟨β⟩ ] Mean of the seasonally forced transmission rate β(t) expressed per unit Δt per susceptible per infected.

alpha

[ α ] Amplitude of the seasonally forced transmission rate β(t) relative to the mean.

epsilon

[ ϵ ] Standard deviation of the random phase shift in the seasonally forced transmission rate β(t).

n

Numeric scalar. Time between the first and last observations in units Δt, so that simulated time series have n+1 observations.

with_dem_stoch

Logical scalar. If TRUE, then the simulation is generated by adaptivetau::ssa.adaptivetau(). Otherwise, it is generated by deSolve::ode() (see Details).

seed

Integer scalar or NA. If an integer, then seed is passed to set.seed() internally, making the simulation, which is randomly generated (see Details), reproducible with an identical call to make_data().

ode_control

A list of arguments to be passed to deSolve::ode(), specifying options for numerical integration, such as method, rtol, and atol. Not used if with_dem_stoch = TRUE.

Value

A data frame with numeric columns:

t

Time in units Δt. Equal to t0 + seq(0, n, by = 1).

t_years

Time in years. Equal to (t0 + seq(0, n, by = 1)) * dt_weeks * (7 / 365).

beta

Seasonally forced transmission rate expressed per unit Δt, per susceptible per infected, without environmental noise.

beta_phi

Seasonally forced transmission rate expressed per unit Δt per susceptible per infected, with environmental noise.

N

Population size.

S

Number of susceptibles.

I

Number of infecteds.

R

Number of removeds.

Q

Cumulative incidence. Q[i] is the number of infections between times t[1] and t[i].

Z

Incidence. Z[i] is the number of infections between times t[i-1] and t[i], computed as Q[i] - Q[i-1].

C

Reported incidence. C[i] is the number of infections reported between times t[i-1] and t[i], equal to the number of successes in Z[i-round(trep)] independent Bernoulli trials, with success probability prep.

A list of the arguments of make_data() is included as an attribute.

Model

make_data() generates data at times tk = t0+kΔt (for k = 0,...,n) using the system of SIR equations

S′ = νcÑ0β(t)SIμcS
I′ = β(t)SIγIμcI
R′ = γIμcR

with γ = 1 / tgen and β(t) = ⟨β⟩ (1 + α cos(2πt / (1 year) + ϕ(t;ϵ))), where ϕ(t;ϵ) is the linear interpolant of random noise {(tk,Φk)} with Φk ~ Normal(0,ϵ2).

Randomness

Simulations have three sources of randomness:

  1. Environmental stochasticity

    • The seasonal forcing function contains a randomly generated phase shift ϕ(t;ϵ). Randomness is introduced by choosing epsilon greater than 0.

  2. Demographic stochasticity

    • If with_dem_stoch = TRUE, then observations are generated by realizing a continuous-time stochastic process, in which event probabilities are determined by terms in the ODE. See adaptivetau::ssa.adaptivetau().

    • Probability of infected decrease is set to zero whenever the number of infecteds is 1, in order to prevent disease fadeout.

    • If with_dem_stoch = FALSE, then observations are generated by numerically integrating the ODE. See deSolve::ode().

  3. Under-reporting of cases

    • Reported incidence C is obtained from incidence Z by modeling C[i] as the number of successes in Z[i-round(trep)] independent Bernoulli trials, with success probability prep. Randomness is introduced by choosing prep in (0,1).

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
# Deterministic simulation
par_list <- make_par_list(dt_weeks = 1)
df <- make_data(
  par_list = par_list,
  n = 20 * 365 / 7,       # 20 years is ~1042 weeks
  with_dem_stoch = FALSE  # numerical integration of d(S,I,R)/dt
)
head(df)

# Reproducible stochastic simulation
par_list <- make_par_list(
  dt_weeks = 1,
  epsilon = 0.5, # s.d. of environmental noise
  prep = 0.5     # case reporting probability
)
df <- make_data(
  par_list = par_list,
  n = 20 * 365 / 7,      # 20 years is ~1042 weeks
  with_dem_stoch = TRUE, # stochastic simulation of d(S,I,R)/dt
  seed = 5               # seed for RNG
)
head(df)

davidearn/fastbeta documentation built on June 14, 2020, 3:11 p.m.