ptpi: Peak to Peak Iteration

ptpiR Documentation

Peak to Peak Iteration

Description

Estimates the initial sizes of the susceptible, infected, and removed populations from a periodic, equally spaced incidence time series and other data. Interpret with care, notably when supplying time series that are only “roughly” periodic.

Usage

ptpi(series, constants, a = 0L, b = nrow(series) - 1L,
     tol = 1e-03, iter.max = 32L,
     complete = FALSE, backcalc = FALSE, ...)

Arguments

series

a “multiple time series” object, inheriting from class mts, with three columns storing (“parallel”, equally spaced) time series of incidence, births, and the per capita natural mortality rate, in that order.

constants

a numeric vector of the form c(Sa, Ia, Ra, gamma, delta), specifying a starting value for the state at time a and rates of removal and loss of immunity, in that order.

a

the time of the first peak in the incidence time series. It is rounded internally to generate a 0-index of rows of series.

b

the time of the last peak in the incidence time series that is in phase with the first. It is rounded internally to generate a 0-index of rows of series.

tol

a tolerance indicating a stopping condition; see ‘Details’.

iter.max

the maximum number of iterations.

complete

a logical indicating if the result should preserve the state at times a, ..., b in each iteration.

backcalc

a logical indicating if the state at time 0 should be back-calculated.

...

optional arguments passed to deconvolve, if the first column of series represents reported incidence or mortality rather than incidence.

Details

ptpi works by computing the iteration described in fastbeta from time t = a to time t = b, iteratively, until the relative difference between the states at times a and b descends below a tolerance. The state at time a in the first iteration is specified by the user. In subsequent iterations, it is the value of the state at time b calculated in the previous iteration.

If backcalc = FALSE, then ptpi returns a list with component value equal to the state at time b in the last iteration. By periodicity, this value estimates the “true” state at time a.

If backcalc = TRUE, then value is back-calculated so that it estimates the “true” state at time 0. This works by inverting the transformation defining one step of the iteration, hence (components of) the back-calculated result can be nonsense if the inverse problem is ill-conditioned.

Value

A list with elements:

value

the estimated value of the initial state, which is the state at time a if backcalc = FALSE and the state at time 0 if backcalc = TRUE.

delta

the relative difference computed in the last iteration.

iter

the number of iterations performed.

X

if complete = TRUE, then a “multiple time series” object, inheriting from class mts, with dimensions c(b-a+1, 3, iter). X[, , i] gives the state at times a, ..., b in iteration i.

References

Jagan, M., deJonge, M. S., Krylova, O., & Earn, D. J. D. (2020). Fast estimation of time-varying infectious disease transmission rates. PLOS Computational Biology, 16(9), Article e1008124, 1-39. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1371/journal.pcbi.1008124")}

Examples


data(sir.e01, package = "fastbeta")
a <- attributes(sir.e01)
str(sir.e01)
plot(sir.e01)

## We suppose that we have perfect knowledge of incidence,
## births, and the data-generating parameters, except for
## the initial state, which we "guess"
series <- cbind(sir.e01[, c("Z", "B")], mu = a[["mu"]](0))
colnames(series) <- c("Z", "B", "mu") # FIXME: stats:::cbind.ts mangles dimnames
constants <- c(Sa = sir.e01[[1L, "S"]],
               Ia = sir.e01[[1L, "I"]],
               Ra = sir.e01[[1L, "R"]],
               gamma = a[["gamma"]],
               delta = a[["delta"]])

plot(series[, "Z"])
a <- 8; b <- 216
abline(v = c(a, b), lty = 2)

L <- ptpi(series, constants, a = a, b = b, complete = TRUE, tol = 1e-06)
str(L)

S <- L[["X"]][, "S", ]
plot(S, plot.type = "single")
lines(sir.e01[, "S"], col = "red", lwd = 4) # the "truth"
abline(h = L[["value"]]["S"], v = a, col = "blue", lwd = 4, lty = 2)

## The relative error
L[["value"]]["S"] / sir.e01[1L + a, "S"] - 1

davidearn/fastbeta documentation built on April 30, 2024, 2:35 a.m.