---
title: Reproductive Number 1 - Practice
output:
  html_document:
    theme: null
    highlight: null
    fig_retina: null
    fig_caption: true
    mathjax: default 
    keep_md: false
bibliography: references.bib
---


```{r, echo = FALSE}
#this code loads the settings file for the current app so we can automatically 
#list the functions in the further information section
currentrmdfile = knitr::current_input() 
currentappinfo = gsub("_documentation.Rmd" ,"_settings.R",currentrmdfile)
source(currentappinfo)
```


## Overview {#shinytab1}
This app is meant to teach you about the basic concepts behind the reproductive number. Read about the model in the "Model" tab. Then do the tasks described in the "What to do" tab.

### Learning Objectives
* Understand the concept of the reproductive number
* Be able to compute the reproductive number for simple models
* Be able to determine the reproductive number using final outbreak size information
* Be able to determine the reproductive number using initial infected numbers


## The Model {#shinytab2}

### Model Overview
For this app, we'll use the basic compartmental SIR model. We allow for 3 different stages/compartments: 

* **S** - uninfected and susceptible individuals 
* **I** - infected and infectious individuals (note that these terms are often used interchangeably, but technically we are talking about someone who is infected **and** is infectious, i.e. can infect others)
* **R** - recovered/removed individuals. Those are individuals that do not further participate, either because they are now immune or because they died.


In addition to specifying the *compartments* of a model, we need to specify the dynamics determining the changes for each compartment. Broadly speaking, some processes increase the number of individuals in a given compartment/stage and other processes that lead to a reduction. Those processes are sometimes called inflows and outflows. 

For our system, we specify the following processes/flows: 

* Susceptible individuals (_S_) can become infected by infectious individuals (_I_) at some rate, _b_. This leads to the susceptible individual leaving the _S_ compartment and entering the _I_ compartment.   
* Infected individuals recover and enter the recovered (_R_) compartment at some rate, _g_.



### Model Implementation
The flow diagram and the set of equations which are used are the basic SIR model and are shown again:



```{r diagram, fig.cap='Flow diagram for simple SIR model.',  echo=FALSE}
knitr::include_graphics("../media/basicSIRmodelfigure.png")
```

Model equations:

$$\dot S = - b SI $$ 
$$\dot I = b S I - g I $$
$$\dot R = g I $$

The underlying model implements one other process. In previous apps, you learned that an unrealistic feature of models implemented by ordinary differential equations (ODE) is that the number of individuals can drop below 1, i.e. there could be fraction of infected. For this app, the underlying code is such that if the number of infected drops below 1, it is set to 0. This is in some way a bit of a "hack" to deal with this issue. When you work through the stochastic apps, you will learn better ways of treating this. The advantage of doing it the way it's done here is that we can keep using the ODE model formulation, without stochasticity, which makes things easy.


### Reproductive number

_The app and tasks deal with the reproductive number concept. The following section provides a very brief introduction. I recommend reading a bit more about it. I'm following the terminology of [my own write-up](https://ahgroup.github.io/DSAIDEbook/reproductive-number.html). You can also check the books listed in the ID introduction app or some of the papers listed in the _Further Resources_ section of this app._

The reproductive number is defined as the average number of new infected (and infectious) individuals caused by one infectious individual. The _basic_ reproductive number is the reproductive number in a scenario where everyone is susceptible. For the SIR model shown above, this quantity is 
$$R_0=\frac{bS_0}{g}$$ 
where _S~0~_ is the initial number of susceptibles.

For a single outbreak (no births, natural death or waning immunity) and a basic SIR model, an equation linking the final number of susceptibles left at the end and the basic reproductive number is 
$$R_0=\frac{\ln(S_f)}{(S_f - 1)}$$
where $\ln()$ is the natural logarithm and _S~f~_ is the __fraction__ of susceptibles still left.


Note the unfortunate fact that the letter _R_ is used both for the recovered compartment in the model and the reproductive number. This is standard notation and I'll therefore use it here. Just be careful to figure out from the context if someone is talking about the recovered individuals or the reproductive number.

## What to do {#shinytab3}


**The tasks below are described in a way that assumes everything is in units of months (rate parameters, therefore, have units of inverse months). If any quantity is not given in those units, you need to convert it first (e.g. if it says a year, you need to convert it to 12 months).**

### Task 1
* Set the simulation with 1000 susceptibles and 1 infected.
* Simulation time 12 months, _g_=5, _b_=0.01. 
* Run the simulation, you should get an outbreak. Use the final size equation linking R~0~ and the fraction of susceptible hosts left at the end of the outbreak to compute the reproductive number (see the information in the _Model_ tab). 

### Task 2
* Use the equation that expresses R~0~ as a function of model parameters for the simple SIR model. Using the values of the model parameters, compute R~0~ that way. Check that it agrees with what you found in the previous task.

### Task 3 
* Double the value of the transmission parameter, _b_. Leave everything else as before.
* Before you run the simulation, what do you expect to see and what do you expect to get for R~0~?
* Run the simulation and compute R~0~ using the final outbreak size to test your expectations. 

### Task 4
* Double the rate of the recovery parameter, _g_. Leave everything else unchanged.
* Think about your expectations for R~0~ and the resulting outbreak dynamics.
* Run the simulation to check your expectations. Use the final outbreak size to compute R~0~.

### Task 5 
* Set the model parameters back to those given in task #1.
* Another way to estimate R~0~ is to determine the rate of increase in infected hosts at the beginning of the outbreak. During the initial phase, new infections increase exponentially according to _I(t)=I~0~ exp(rt)_, with _r_ being the rate of growth. Usually, for any real outbreak, you do not know the number of infected at the start, I~0~, or the exact time the outbreak starts. It is still possible to estimate _r_ by obtaining two values of _I_ at two time points during that initial growth rate, i.e. _I~1~_ at time _t~1~_ and _I~2~_ at time _t~2~_. One obtains equation _I~1~=I~0~ exp(r t~1~)_ and _I~2~=I~0~ exp(r t~2~)_. By solving one of these equations for _I~0~_ and substituting into the other, we get _I~2~= I~1~ exp(r (t~2~ - t~1~))_. By solving the model for _r_ and entering numbers for _I~1~_ and _I~2~_ and times _t~1~_ and _t~2~_ we can figure out _r_.
* Let's try that. Run the model with tmax = 0.2 and tmax = 0.4 and record the number of infected at the end of the simulation for each time. Then substitute all the values into the equation you found for _r_ and thus compute the growth rate.
* For this model, the growth rate and R~0~ are related through _R~0~ = 1+rD_, where _D_ is the average duration of the infectious period (i.e. the inverse of the recovery rate). Use this to determine R~0~. You should get essentially the same answer (up to some rounding differences) as for task #1.
* Note that the choice of _t~1~_ and _t~2~_ can influence the results. Earlier times are better since once the number of susceptibles starts to drop markedly, the growth of infected slows down and is not exponential anymore.


### Task 6 
* What is the value of the reproductive number _R_ at the time the outbreak peaks? (It's only called R~0~ at the beginning for a fully susceptible population). Explain how you can find that value for R, both using intuitive reasoning and using the equation for R~0~ given above (R~0~ = 1+rD). Note that at this R value, the outbreak wanes, but people still get infected. What R value would you need to halt any further infections completely?

### Task 7
* What would happen if a new ID came along that had an R~0~ value that was the same as the one you just determined in the previous question, namely the value of R at the peak of an outbreak? Test this with the simulation. Choose parameter values for which you get the right R~0~ and observe what the model produces.


### Task 8 
* R~0~ quantifies the level of transmissibility of an ID, which determines how many people will become infected or what level of intervention is needed to stop/prevent an outbreak. However, it is important to be aware that R~0~ says nothing about the timing/dynamics of the outbreak. Set parameter values as in #1. Run an outbreak, record the time of peak and duration of the outbreak (the latter is somewhat ill-defined, so just come up with a rough number). Then increase the infectious duration by a factor of 4 (rate reduced by a factor of 4) and adjust the infectiousness-level such that you get the same R~0~ as before. Run again and compare the results concerning total outbreak size and timing of outbreak.    






## Further Information {#shinytab4}

* This app (and all others) are structured such that the Shiny part (the graphical interface you see and the server-side function that goes with it) calls an underlying R script (or several) which runs the simulation for the model of interest and returns the results.
* For this app, the underlying function running the simulation is called ``r appsettings$simfunction``. You can call them directly, without going through the shiny app. Use the `help()` command for more information on how to use the functions directly. If you go that route, you need to use the results returned from this function and produce useful output (such as a plot) yourself. 
* You can also download all simulator functions and modify them for your own purposes.  Of course to modify these functions, you'll need to do some coding.
* For examples on using the simulators directly and how to modify them, read the package vignette by typing `vignette('DSAIDE')` into the R console.
* For the equations mentioned in the tasks, see e.g. [@heffernan05a; @keeling08; @vynnycky10].
* For some information on the reproductive number, see e.g. [@fine11; @heffernan05a; @roberts07; @wallinga07].

### References



