This html notebook describes the COVID-19 breakthrough case surveillance project and integrates technical definitions, underlying R code, and corresponding outputs. Navigate to sections of interest by clicking on the table of contents on the left. You can choose to show or hide the code by toggling the code button near the top right of relevant sections. This notebook is divided into four sections, namely, “Introduction”, “Technical notes”, “Data prep”, and “Outputs”.

The contents of this notebook are for internal BCCDC use only. For the most up-to-date public-facing data on COVID-19 breakthrough cases, please refer to the BCCDC Regional Surveillance Dashboard.

Content

Three data sources are used in the surveillance of breakthrough cases:

• BCCDC Public Health Reporting Data Warehouse (PHRDW) integrated COVID dataset, which includes the Health Authority case linelist. These data are used to primarily to assess counts and numerator in rates.

• BCCDC combined population file, which houses BC Stats Population Estimates for previous years, and PEOPLE population projection for current and future years. These data, in combination wih the PHSA Provincial Immunization Registry, are used to assign a vaccination status to people in BC.

Content

• Rates are presented as daily incidence proportion, in which the numerator is the number of new daily cases in each vaccination status group, and the denominator is the number of people in each vaccination status group for each day
• As the immunization program progresses, the number of people in each vaccination status group changes every day
• The PHSA Provincial Immunization Registry and the BCCDC combined population file are used to estimate the number of people in each group in each day
  • Vaccinated, 3 doses = cumulative number of 3rd doses administered
  • Vaccinated, 2 doses = cumulative number of 2nd doses administered - Vaccinated, 3 doses
  • Vaccinated, 1 doses = cumulative number of 1st doses administered - (Vaccinated, 2 doses + Vaccinated, 3 doses)
  • Unvaccinated = population denominator - (Vaccinated, 1 dose + Vaccinated, 2 doses + Vaccinated, 3 doses)
• The dates are shifted to construct a “protected date” concept, in which the breakthrough case numerators and denominators for rates are aligned. The protected date is calculated for each group as follows:
  • “Vaccinated, 1 dose”: 1st dose administration date + 21 days
  • “Vaccinated, 2 doses”: 2nd dose administration date + 14 days
  • “Vaccinated, 3 doses”: 3rd dose administration date + 14 days
• The BC population estimates change from year to year, and this is most noticeable in the “Unvaccinated” strata in the figure below (e.g. See increases and decreases starting 2022)

vaccine_denom_bc <- vaccine_data_clean %>% 
  group_by(protected_date_14, series_dose_number_derived) %>%
  summarise(vaccine_denom = n()) %>%
  ungroup() %>%
  group_by(series_dose_number_derived) %>%
  pad() %>%
  fill_by_value(vaccine_denom, 0) %>%
  ungroup() %>%
  mutate(vax_flag = case_when(series_dose_number_derived == 1 ~ "Vaccinated, 1 dose", 
                              series_dose_number_derived == 2 ~ "Vaccinated, 2 doses", 
                              series_dose_number_derived == 3 ~ "Vaccinated, 3 doses")) %>%
  select(-series_dose_number_derived) %>%
  complete(protected_date_14, vax_flag, fill = list(vaccine_denom = 0)) %>%
  group_by(vax_flag) %>%
  arrange(protected_date_14) %>%
  mutate(cum_vaccine_denom = cumsum(vaccine_denom)) %>%
  ungroup() %>%
  select(-vaccine_denom) %>%
  spread(vax_flag, cum_vaccine_denom) %>%
  mutate(year = year(protected_date_14)) %>%
  left_join(bc_pop_denom_overall) %>%
  mutate(`Vaccinated, 2 doses` = `Vaccinated, 2 doses` - `Vaccinated, 3 doses`,
         `Vaccinated, 1 dose` = `Vaccinated, 1 dose` - (`Vaccinated, 2 doses` + `Vaccinated, 3 doses`), 
         Unvaccinated = pop_denom - (`Vaccinated, 1 dose` + `Vaccinated, 2 doses` + `Vaccinated, 3 doses`)) %>%
  select(-c(pop_denom, year)) %>%
  gather("vax_flag", "cum_vaccine_denom", 2:5)

denom_age_plot_data <- cases_hosp_breakthrough_master_data %>%
  filter(date >= ymd("2021-07-01") & 
           surveillance_ha == "BC" &
           !age_group == "all ages" &
           !is.na(age_group) & 
           case_hosp_death == "cases" &
           metric == "cum_vaccine_denom") %>% 
  arrange(match(age_group, age_group_labels)) %>%
  mutate(age_group = as_factor(age_group))

denom_age_plot <- denom_age_plot_data %>%
  ggplot(aes(x = date, y = value, group = vax_flag_14, color = vax_flag_14)) +
  geom_line(size = 1, 
            show.legend = FALSE) +
  geom_point(data = denom_age_plot_data %>%
               group_by(vax_flag_14, age_group) %>%
               filter(date == max(date)) %>%
               ungroup(), 
             size = 3, 
             show.legend = TRUE) +
  facet_wrap(~age_group) +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 0, vjust = 0.5, hjust = 0.5, size = 15), 
        axis.text.y = element_text(size = 20),
        axis.title.x = element_blank(),  
        axis.title.y = element_text(size = 20, margin = margin(t = 0, r = 20, b = 0, l = 0)), 
        plot.caption = element_text(size = 11, margin = margin(t = 20, r = 0, b = 0, l = 0), hjust = 0),
        panel.grid.minor = element_blank(), 
        panel.grid.minor.x = element_blank(),
        panel.grid.major.x = element_blank(),
        panel.spacing = unit(2, "lines"),
        strip.text = element_text(size = 20), 
        legend.position = "top", 
        legend.title = element_blank(), 
        legend.text = element_text(size = 15, margin = margin(r = 10, unit = "pt")), 
        axis.ticks = element_line()) +
  ylab("Number of people") +
  scale_color_manual(values = plot_colors) +
  scale_y_continuous(label = unit_format(scale = 0.001, unit = "K")) +
  labs(caption = paste0("Data extracted from PHSA Provincial Immunization Registry ", format(today, "%d %b %Y"))) +
  scale_x_date(breaks = trend_date_breaks, labels = trend_date_labels, expand = c(0.05, 0.05))

denom_age_plot