rtestim

https://github.com/dajmcdon/rtestim

Repository

README.Rmd

---
output: github_document
---



```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = FALSE,
  comment = "#>",
  fig.path = "man/figures/README-",
  fig.align = "center",
  fig.height = 3,
  fig.width = 6,
  dpi = 300,
  out.width = "80%"
)
```

# rtestim 


[![R-CMD-check](https://github.com/dajmcdon/rtestim/actions/workflows/R-CMD-check.yaml/badge.svg)](https://github.com/dajmcdon/rtestim/actions/workflows/R-CMD-check.yaml)
[![CRAN status](https://www.r-pkg.org/badges/version/rtestim)](https://CRAN.R-project.org/package=rtestim)


This package uses Poisson likelihood with a trend filtering penalty (a type of 
regularized nonparametric regression)
to estimate the effective reproductive number, Rt.
This value roughly says "how many new infections will result from
each new infection today". Values larger than 1 indicate that an 
epidemic is growing while those less than 1 indicate decline.

## Installation

You can install the development version of rtestim from 
[GitHub](https://github.com/) with:

``` r
# install.packages("remotes")
remotes::install_github("dajmcdon/rtestim")
```

Or the released version on CRAN

``` r
install.packages("rtestim")
```

## Quick example

Here we create some data that "looks" like a typical wave in an epidemic.
Because the model uses regularized regression, we estimate the model
at a range of tuning parameters simultaneously. 

```{r plot-data, message=FALSE, fig.align='center'}
set.seed(12345)
library(rtestim)
library(ggplot2)
dat <- data.frame(
  Time = 1:101,
  incident_cases = c(1, rpois(100, dnorm(1:100, 50, 15) * 500 + 1))
)
ggplot(dat, aes(Time, incident_cases)) +
  geom_point(colour = "cornflowerblue") +
  theme_bw()
```

We fit the model and visualize the resulting estimated sequences of $R_t$:

```{r full-fit, fig.align='center'}
mod <- estimate_rt(observed_counts = dat$incident_cases, nsol = 20)
plot(mod)
```

The additional parameter `nsol = 20` specifies the number of 
tuning parameters for which $R_t$ is estimated. A built in function for 
cross-validation can be used to select the tuning parameter.

```{r cv-estimate, fig.align='center'}
mod_cv <- cv_estimate_rt(dat$incident_cases, nsol = 20)
plot(mod_cv, which_lambda = "lambda.1se")
```

Owner

GitHub Events

Packages