estimRmcmc {EpiLPS}  R Documentation 
Estimation of the reproduction number with LaplacianPsplines via MCMC
Description
This routine estimates the instantaneous reproduction number R_t
;
the mean number of secondary infections generated by an infected individual
at time t
(White et al. 2020); by using Bayesian Psplines and Laplace
approximations (Gressani et al. 2022). The inference approach is fully
stochastic with a Metropolisadjusted Langevin algorithm. The
estimRmcmc()
routine estimates R_t
based on a time series of
incidence counts and a (discretized) serial interval distribution. The
negative binomial distribution is used to model incidence count data and
Psplines (Eilers and Marx, 1996) are used to smooth the epidemic curve.
The link between the epidemic curve and the reproduction number is
established via the renewal equation.
Usage
estimRmcmc(incidence, si, K = 30, dates = NULL, niter = 5000, burnin = 2000,
CoriR = FALSE, WTR = FALSE, priors = Rmodelpriors(), progressbar = TRUE)
Arguments
incidence 
A vector containing the incidence time series. If

si 
The (discrete) serial interval distribution. 
K 
Number of Bsplines in the basis. 
dates 
A vector of dates in format "YYYYMMDD" (optional). 
niter 
The number of MCMC samples. 
burnin 
The burnin size. 
CoriR 
Should the 
WTR 
Should the 
priors 
A list containing the prior specification of the model hyperparameters as set in Rmodelpriors. See ?Rmodelpriors. 
progressbar 
Should a progression bar indicating status of MCMC algorithm be shown? Default is TRUE. 
Value
A list with the following components:
incidence: The incidence time series.
si: The serial interval distribution.
RLPS: A data frame containing estimates of the reproduction number obtained with the LaplacianPsplines methodology.
thetahat: The estimated vector of Bspline coefficients.
Sighat: The estimated variancecovariance matrix of the Laplace approximation to the conditional posterior distribution of the Bspline coefficients.
RCori: A data frame containing the estimates of the reproduction obtained with the method of Cori (2013).
RWT: A data frame containing the estimates of the reproduction obtained with the method of WallingaTeunis (2004).
LPS_elapsed: The routine real elapsed time (in seconds) when estimation of the reproduction number is carried out with LaplacianPsplines.
penparam: The estimated penalty parameter related to the Pspline model.
K: The number of Bsplines used in the basis.
NegBinoverdisp: The estimated overdispersion parameter of the negative binomial distribution for the incidence time series.
optimconverged: Indicates whether the algorithm to maximize the posterior distribution of the hyperparameters has converged.
method: The method to estimate the reproduction number with LaplacianPsplines.
optim_method: The chosen method to to maximize the posterior distribution of the hyperparameters.
HPD90_Rt: The
90\%
HPD interval for Rt obtained with the LPS methodology.HPD95_Rt: The
95\%
HPD interval for Rt obtained with the LPS methodology.
Author(s)
Oswaldo Gressani oswaldo_gressani@hotmail.fr
References
Gressani, O., Wallinga, J., Althaus, C. L., Hens, N. and Faes, C. (2022). EpiLPS: A fast and flexible Bayesian tool for estimation of the timevarying reproduction number. Plos Computational Biology, 18(10): e1010618.
Cori, A., Ferguson, N.M., Fraser, C., Cauchemez, S. (2013). A new framework and software to estimate timevarying reproduction numbers during epidemics. American Journal of Epidemiology, 178(9):1505–1512.
Wallinga, J., & Teunis, P. (2004). Different epidemic curves for severe acute respiratory syndrome reveal similar impacts of control measures. American Journal of Epidemiology, 160(6), 509516.
White, L.F., Moser, C.B., Thompson, R.N., Pagano, M. (2021). Statistical estimation of the reproductive number from case notification data. American Journal of Epidemiology, 190(4):611620.
Eilers, P.H.C. and Marx, B.D. (1996). Flexible smoothing with Bsplines and penalties. Statistical Science, 11(2):89121.
Examples
# Illustration on the 2009 influenza pandemic in Pennsylvania.
data(influenza2009)
epifit_flu < estimRmcmc(incidence = influenza2009$incidence, dates = influenza2009$dates,
si = influenza2009$si[1], niter = 2500,
burnin = 1500, progressbar = FALSE)
tail(epifit_flu$RLPS)
summary(epifit_flu)
plot(epifit_flu)