```{r setup, include=FALSE}
knitr::opts_chunk$set(
  echo = FALSE,
  eval = FALSE,
  cache = FALSE,
  warning = FALSE,
  message = FALSE
)
options(width = 80, digits = 3)
```  
  
# `r emo::ji("target")` Objectives

In this tutorial, you will learn

- load in and manipulate COVID-19 data
- write a model to estimate the growth rate of PCR-detected cases
- compare growth rates in LGAs with vaccination coverage

# `r emo::ji("wrench")` Preparation

Data links:

- [Case data](https://gist.githubusercontent.com/MikeLydeamore/e911c95a792c1aeb72368a044de5ed27/raw/c8c9fe8842b25353d0ae60d447d8cb0a157604af/tutorial_10_vic_cases_by_lga.csv)

- [Vaccination data](https://gist.githubusercontent.com/MikeLydeamore/d445691db1922b355569e45824d46f42/raw/3bc674501c37f412fa2f32711d60f658b01638ca/tutorial_10_vic_vaccinations_by_lga.csv)

To load in the data, we can read it directly from the URLs above.
```{r, message=F}
case_data <- readr::read_csv("https://gist.githubusercontent.com/MikeLydeamore/e911c95a792c1aeb72368a044de5ed27/raw/c8c9fe8842b25353d0ae60d447d8cb0a157604af/tutorial_10_vic_cases_by_lga.csv")

vaccination_data <- readr::read_csv("https://gist.githubusercontent.com/MikeLydeamore/d445691db1922b355569e45824d46f42/raw/3bc674501c37f412fa2f32711d60f658b01638ca/tutorial_10_vic_vaccinations_by_lga.csv")
```

We'll also need the `tidyverse` and `lubridate`:
```{r load-packages}
library(tidyverse)
library(lubridate)
```
# `r emo::ji("gear")` Tutorial

## 1. Manipulating the data

a. Convert the long `case_data` data frame into a table with three columns:

    i) `date`: The date of diagnosis
    ii) `local_government_area`: The LGA of interest
    iii) `num_cases`: The number of PCR positive cases on a given day.

```{r aggregate-cases, message=F}
case_data_aggregated <- case_data %>%
  group_by(diagnosis_date, Localgovernmentarea) %>%
  summarise(num_cases = n()) %>%
  rename(
    date = diagnosis_date,
    local_government_area = Localgovernmentarea
  )

head(case_data_aggregated)
```


b. Draw a plot with the _total_ number of cases on the y axis, and the LGAs (ordered from least to most cases) on the x axis.  
You can add `theme(axis.text.x = element_text(angle = 90, hjust=1, vjust=0.5))` to your `ggplot` to rotate and align the text.

c. Convert the dose numbers for vaccination to numeric variables instead of characters. For LGAs with >95% coverage, set them to 0.95.  
(hint: `substr(x, 1, nchar(x)-1)` will remove the last character of a string)

```{r convert-vax-to-numeric, message=F}
vaccination_data_converted <- vaccination_data %>%
  mutate(
    FIRST = substr(FIRST, 1, nchar(FIRST) - 1),
    SECOND = substr(SECOND, 1, nchar(SECOND) - 1),
    THIRD = substr(THIRD, 1, nchar(THIRD) - 1),
  ) %>%
  mutate(
    FIRST = case_when(
      FIRST == ">95" ~ 0.95,
      TRUE ~ as.numeric(FIRST) / 100
    ),
    SECOND = case_when(
      SECOND == ">95" ~ 0.95,
      TRUE ~ as.numeric(SECOND) / 100
    ),
    THIRD = case_when(
      THIRD == ">95" ~ 0.95,
      TRUE ~ as.numeric(THIRD) / 100
    )
  ) %>%
  rename(
    local_government_area = LGA
  )

head(vaccination_data_converted)
```

d. Draw a plot with the third dose coverage on the y axis, and LGAs (ordered from lowest to highest coverage) on the x axis.

## 2. Calculating the growth rate

The growth rate over a time period, $t_1$, is defined as

$$C(t_1) = C(t_0) \exp(rt),$$ where $r$ is the rate of growth.

Rearranging this expression gives,
$$r = \frac{\log\left(\frac{C(t_1)}{C(t_0)}\right)}{t}.$$

Below is a function called `calculate_growth_rate` that takes three inputs (`Ct0`, `Ct1`, and `time`) and returns the growth rate.  
Remember, it's always a good idea to test your code and make sure it returns the numbers you expect.

```{r growth-rate-function, echo=TRUE}
calculate_growth_rate <- function(Ct0, Ct1, time) {
  log(Ct1 / Ct0) / time
}
```

a. Calculate the 7-day growth rate for each Victorian LGA since 1 January 2022. For this tutorial, we will assume detection date (which is what is present in the data) is the infection date to calculate the growth rate on.    

```{r calculate-growth-rate, message=F}
growth_rate <- case_data_aggregated %>%
  filter(year(date) >= 2022) %>%
  group_by(local_government_area) %>%
  mutate(
    growth_rate = slider::slide_dbl(
      .x = num_cases,
      .f = function(v) {
        calculate_growth_rate(Ct0 = v[1], Ct1 = v[7], time = 7)
      },
      .before = 6,
      .after = 0
    )
  )

head(growth_rate)
```

## 3. Correlating growth rates with vaccination coverage

You've been asked to explore whether vaccination coverage and growth rate are correlated. You are on a tight deadline and so only have time to produce one plot.

Thankfully, your colleague has given you some code to join together the vaccination coverage and growth rate data:

```{r remove-suffixes, echo=TRUE}
codes_to_remove <- c(
  " \\(C\\)",
  " \\(S\\)",
  " \\(RC\\)",
  " \\(B\\)"
) %>%
  paste0(collapse = "|")

summarised_growths <- summarised_growths %>%
  mutate(local_government_area = gsub(
    pattern = codes_to_remove,
    replacement = "",
    x = local_government_area
  ))
```

```{r third-dose-join, echo=TRUE}
terminal_coverage <- vaccination_data_converted %>%
  left_join(summarised_growths, on = "local_government_area")
```

Create a plot that shows evidence _for_ or _against_ the hypothesis that vaccination coverage and growth rate are correlated. Write a maximum four sentence executive summary of your findings, and compare with a partner in the class.

#### © Copyright 2022 Monash University