diff --git a/maternal_VE_study/BronchSTOP_analysis_for_submission.qmd b/maternal_VE_study/BronchSTOP_analysis_for_submission.qmd
deleted file mode 100644
index 9d8b5d455de211dd38638f75a02f7f5af571687e..0000000000000000000000000000000000000000
--- a/maternal_VE_study/BronchSTOP_analysis_for_submission.qmd
+++ /dev/null
@@ -1,769 +0,0 @@
----
-title: "BronchStop data quality report and analysis"
-format: pdf
----
-
-## Analysis on `r format(Sys.Date(),'%d/%m/%Y')`
-
-```{r setup}
-#| echo: false
-#| message: false
-
-library(tidyverse)
-library(janitor)
-library(tidyr)
-library(glue)
-library(gt)
-library(survival)
-library(broom)
-library(gtsummary)
-
-load(file = "bronchstop_livedata.rdata")
-
-print(glue("Number of consented subjects in database: {data%>%filter(consented==TRUE)%>%nrow()}"))
-
-data<-data%>%
- filter(consented==TRUE)%>%
- mutate(vax_gest=case_when(vax_gest %in% c(1:18) ~ vax_gest+23,
- vax_gest == 99 ~ 99))%>%
- mutate(vax_gest_record=case_when(vax_gest_record %in% c(1:18) ~ vax_gest_record+23,
- vax_gest_record == 99 ~ 99))%>%
- mutate(age_gest=case_when(age_gest %in% c(1:15) ~ age_gest+21,
- age_gest == 99 ~ 99))%>%
- mutate(age_gest=replace_na(age_gest,40))%>%
- mutate(age_gest_term_or_not=if_else(age_gest<36,FALSE,TRUE))%>%
- mutate(vax_gest_born=age_gest-vax_gest)%>%
- mutate(potential_excludes=if_else(vax_gest_born %in% c(-1,0,1,2),TRUE,FALSE))%>%
- mutate(LOA = time_length(interval(ymd(attend_dat),ymd(dis_dat_final)), "day"))%>%
- mutate(follow_up_form_complete=case_when(follow_up_form_complete==0 ~"incomplete",
- follow_up_form_complete==1 ~"unverified",
- follow_up_form_complete==2 ~"complete"))%>%
- mutate(vax_yn_record=case_when(vax_yn_record==0 ~"confirmed_unvax",
- vax_yn_record==1 ~"confirmed_vax",
- vax_yn_record==99 ~"checked_not_known",
- is.na(vax_yn_record) ~"not_done_yet"))%>%
- mutate(mum_thinks_vaccinated=case_when(vax_yn==0 ~"thinks_unvax",
- vax_yn==1 ~"thinks_vax",
- vax_yn==99 ~"mum_unsure"))%>%
- mutate(vax_status_combined=case_when(vax_yn_record=="confirmed_vax" ~ TRUE,
- vax_yn_record=="confirmed_unvax" ~ FALSE,
- mum_thinks_vaccinated=="thinks_vax" ~ TRUE,
- mum_thinks_vaccinated=="thinks_unvax" ~ FALSE,
- .default = NA))%>%
- mutate(vax_status=case_when(vax_yn==1 ~ TRUE,
- vax_yn==0 ~ FALSE))%>%
- mutate(pertussis_status=case_when(pert_record==1 ~ TRUE,
- pert_record==0 ~ FALSE))%>%
- #mutate(sex=case_when(sex==2 ~"female",
- # sex==1 ~"male"))%>%
- mutate(male_gender=case_when(sex==2 ~FALSE,sex==1 ~TRUE))%>%
- mutate(age_days = time_length(interval(ymd(dob_infant),ymd(attend_dat)), "day"))%>%
- mutate(age_mnths_calc = time_length(interval(ymd(dob_infant),ymd(attend_dat)), "month"))%>%
- mutate(month=month(attend_dat))%>%
- mutate(age_under_3_months = if_else(age_mnths_calc<4,TRUE,FALSE))%>%
- mutate(case=case_when(virus %in% c("rsv") ~ TRUE,
- virus %in% c("none","other") ~ FALSE))%>%
- rename("comorb-other"=comorb_othr)%>%
- mutate(across(starts_with("comorb_"), ~ as.logical(.x)))%>%
- mutate(across(starts_with("treat_hosp_"), ~ as.logical(.x)))%>%
- mutate(imd_eng=case_when(imd_eng %in% 1:2 ~ 1,
- imd_eng %in% 3:4 ~ 2,
- imd_eng %in% 5:6 ~ 3,
- imd_eng %in% 7:8 ~ 4,
- imd_eng %in% 9:10 ~ 5,
- .default = NA))%>%
- mutate(imd_scot_2=if_else(imd_scot_2==99,NA,imd_scot_2))%>%
- mutate(imd=coalesce(imd_eng,imd_scot_2))%>%
- mutate(breastfed=if_else(breast==99,NA,as.logical(breast)))%>%
- mutate(maternal_ethnicity=case_when(mat_eth==1 ~ "White",
- mat_eth==2 ~ "Asian",
- mat_eth==3 ~ "Black",
- mat_eth==4 ~ "Mixed",
- mat_eth==99 ~ "Other",
- mat_eth==999 ~ NA,
- .default = NA))%>%
- mutate(maternal_white=case_when(mat_eth==1 ~ TRUE,
- mat_eth %in% c(2:99) ~ FALSE,
- mat_eth==999 ~ NA,
- .default = NA))%>%
- mutate(prem=comorb_1)%>%
- filter(attend_dat<dmy("21-01-2025"))
-
-```
-## Data quality check
-### incomplete forms / data errors
-```{r incomplete}
-#| echo: false
-#| message: false
-
-ineligible_exclusions=c("gla-1203","gla-1188","add-580")
-
-#exclude GLA-1203 and GLA-1188
-#exclude add-580 (not admitted)
-cat(glue("Excluding {length(ineligible_exclusions)} who weren't eligible to recruit. ({paste(ineligible_exclusions, collapse=', ')})"))
-data<-data%>%filter(!record_id %in% ineligible_exclusions)
-
-print(glue("Number of remaining subjects: {data%>%filter(consented==TRUE)%>%nrow()}"))
-
-```
-
-\newpage
-### assessing maternal vaccine recall accuracy
-```{r mat_vax_recall_accuracy}
-#| echo: false
-#| message: false
-
-data%>%
- tabyl(vax_yn_record,mum_thinks_vaccinated)%>%
- adorn_totals(where=c("row", "col"))%>%gt()%>%
- tab_header(title = md("validation of maternal recall"))
-
-library(epiR)
-vax_stats<-data%>%
- mutate(pred = case_when(mum_thinks_vaccinated=="thinks_vax" ~ "vax",
- mum_thinks_vaccinated=="thinks_unvax" ~ "unvax"),
- truth = case_when(vax_yn_record=="confirmed_vax" ~ "vax",
- vax_yn_record=="confirmed_unvax" ~ "unvax"))%>%
- tabyl(pred,truth,show_na = F)%>%as.data.frame()%>%arrange(desc(pred))%>%select(vax,unvax)%>%
- as.matrix()
-
-colnames(vax_stats)<-c("recall_vax","recall_unvax")
-rownames(vax_stats)<-c("conf_vax","conf_unvax")
-
-accuracy<-binom.test((vax_stats[1,1]+vax_stats[2,2]),sum(vax_stats))%>%
- tidy()%>%
- select(estimate,conf.low, conf.high)%>%
- mutate(across(everything(),~sprintf("%0.1f%%", .x * 100)))
-
-cat(glue("Accuracy of maternal recall {accuracy$estimate} ({accuracy$conf.low}-{accuracy$conf.high})"))
-
-vax_stats%>%
- epi.tests(., conf.level = 0.95)
-
-```
-
-\newpage
-### Exclusions data
-```{r exclusions}
-#| echo: false
-#| message: false
-
-exclusionTable<-data.frame(reason=NULL,count=NULL)
-# Had monoclonals:
-hadMabs<-data%>%filter(paliz==1 |
- mabs_othr_1==1 |
- mabs_othr_1==2)%>%
- select(record_id)
-cat(glue("Excluding {nrow(hadMabs)} who had had monoclonal antibodies in last month"))
-if(nrow(hadMabs)>0){
- data<-data%>%filter(!record_id %in% hadMabs$record_id)
- exclusionTable<-rbind(exclusionTable, data.frame(reason="had mabs",count=nrow(hadMabs)))
-}
-
-# records not marked as complete:
-notComplete<-data%>%filter(follow_up_form_complete=="incomplete")%>%
- select(record_id)
-
-cat(glue("Excluding {nrow(notComplete)} whos records were not marked as complete"))
-if(nrow(notComplete)>0){
- usefuls<-data%>%filter(record_id %in% notComplete$record_id)%>%select(record_id,case,vax_yn)%>%drop_na(case)%>%filter(!vax_yn==99)
- cat(glue("\t {nrow(usefuls)} of these have enough to be useful though. ({paste(usefuls$record_id, collapse=', ')})"))
- data<-data%>%filter(!record_id %in% notComplete$record_id)
- exclusionTable<-rbind(exclusionTable, data.frame(reason="not completed",count=nrow(notComplete)))
-}
-
-# incomplete testing data:
-incompleteTests<-data%>%filter(!virus %in% c("rsv","none","other"))%>%select(record_id)
-cat(glue("Excluding {nrow(incompleteTests)} who had had no PCR result"))
-if(nrow(incompleteTests)>0){
- data<-data%>%filter(!record_id %in% incompleteTests$record_id)
- exclusionTable<-rbind(exclusionTable, data.frame(reason="not tested",count=nrow(incompleteTests)))
-}
-
-# # mum did not know if vaccinated or not:
- unknownVax<-data%>%filter(is.na(vax_status))%>%select(record_id)
- cat(glue("Excluding {nrow(unknownVax)} with unknown vaccine status"))
- if(nrow(unknownVax)>0){
- data<-data%>%filter(!record_id %in% unknownVax$record_id)
- exclusionTable<-rbind(exclusionTable, data.frame(reason="unknown vax status",count=nrow(unknownVax)))
- }
-
-print(glue("Number of remaining subjects: {data%>%filter(consented==TRUE)%>%nrow()}"))
-
-exclusionTable%>%
- gt()%>%tab_header(title = md("Exclusions"))
-
-data%>%select(record_id,age_mnths_calc,dob_infant,attend_dat)%>%arrange(age_mnths_calc)%>%top_n(5,age_mnths_calc)%>%
- gt()%>%tab_header("oldest kids and their DoBs")
-
-cat("min & max recruitment dates....")
-summary(data$attend_dat)
-
-```
-
-```{r looking_at_gestations}
-#| echo: false
-#| message: false
-#| include: true
-
-#if use recall data as priority
-under2weeks_recall<-
-data%>%
- #calculate differences with gestations
- mutate(vax_gest=if_else(vax_gest==99,NA,vax_gest))%>%
- mutate(vax_gest_record=if_else(vax_gest_record==99,NA,vax_gest_record))%>%
- mutate(age_gest=if_else(age_gest==99,NA,age_gest))%>%
- mutate(vax_gest_combined=coalesce(vax_gest,vax_gest_record))%>%
- mutate(gestVax_to_gestBorn=age_gest-vax_gest_combined)%>%
-
- #calculate differences with dates
- mutate(vax_dat_combined=coalesce(vax_dat,vax_dat_record))%>%
- mutate(dob_to_vaxdate=time_length(interval(ymd(vax_dat_combined),ymd(dob_infant)), "weeks"))%>%
-
- filter(vax_yn==1)%>%
- mutate(dob_to_vaxdate=if_else(dob_to_vaxdate<0,NA,dob_to_vaxdate))%>%
- mutate(vax_diff=coalesce(dob_to_vaxdate,gestVax_to_gestBorn))%>%
- mutate(under2weeks=if_else(vax_diff<2,TRUE,FALSE))%>%
- select(record_id,under2weeks,vax_diff,dob_to_vaxdate,gestVax_to_gestBorn,vax_gest,vax_gest_record,age_gest,vax_dat,vax_dat_record,dob_infant)
-
-
-#if use vax data as priority
-under2weeks_records<-
- data%>%
- #calculate differences with gestations
- mutate(vax_gest=if_else(vax_gest==99,NA,vax_gest))%>%
- mutate(vax_gest_record=if_else(vax_gest_record==99,NA,vax_gest_record))%>%
- mutate(age_gest=if_else(age_gest==99,NA,age_gest))%>%
- mutate(vax_gest_combined=coalesce(vax_gest_record,vax_gest))%>%
- mutate(gestVax_to_gestBorn=age_gest-vax_gest_combined)%>%
-
- #calculate differences with dates
- mutate(vax_dat_combined=coalesce(vax_dat_record,vax_dat))%>%
- mutate(dob_to_vaxdate=time_length(interval(ymd(vax_dat_combined),ymd(dob_infant)), "weeks"))%>%
-
- filter(vax_yn==1)%>%
- mutate(dob_to_vaxdate=if_else(dob_to_vaxdate<0,NA,dob_to_vaxdate))%>%
- mutate(vax_diff=coalesce(dob_to_vaxdate,gestVax_to_gestBorn))%>%
- mutate(under2weeks=if_else(vax_diff<2,TRUE,FALSE))%>%
- select(record_id,under2weeks,vax_diff,dob_to_vaxdate,gestVax_to_gestBorn,vax_gest,vax_gest_record,age_gest,vax_dat,vax_dat_record,dob_infant)
-
-
-under2weeks_recall%>%tabyl(under2weeks)%>%gt()%>%tab_header("prioritising recall to identify <2 weeks since vaccine")
-
-under2weeks_records%>%tabyl(under2weeks)%>%gt()%>%tab_header("prioritising records to identify <2 weeks since vaccine")
-
-data<-data%>%
- left_join(under2weeks_records%>%select(record_id,under2weeks))
-
-```
-
-```{r plotting function}
-#| echo: false
-#| message: false
-#| include: false
-
-calculate_VE<-function(vaxData){
-
- t1<-vaxData%>%
- mutate(case=if_else(case,"case","control"))%>%
- mutate(vax_status=if_else(vax_status,"vax","unvax"))%>%
- tabyl(vax_status,case)%>%
- adorn_totals(where=c("row", "col"))%>%
- rename(" "=vax_status)%>%
- gt()%>%tab_header(title = md("2x2 data feeding into regression"))
-
-tableOutput<-data.frame(formula=NULL,VE=NULL,additional=NULL)
-
-VE_calc<-clogit(formula=case ~ vax_status, data=vaxData)%>%
- tidy(conf.int=TRUE, exponentiate=TRUE)%>%
- mutate(estimate=1-estimate,
- conf_high=1-conf.low,
- conf_low=1-conf.high)%>%
- select(estimate,conf_low,conf_high)
-
-tableOutput<-rbind(tableOutput,data.frame(
- formula="Unadjusted",
- VE=glue("{100 * round(VE_calc$estimate,3)}% ({100 * round(VE_calc$conf_low,3)}% - {100 * round(VE_calc$conf_high,3)}%)"),
- additional="")
-)
-
-VE_calc<-clogit(formula=case ~ vax_status + strata(month), data=vaxData)%>%
- tidy(conf.int=TRUE, exponentiate=TRUE)%>%
- mutate(estimate=1-estimate,
- conf_high=1-conf.low,
- conf_low=1-conf.high)%>%
- select(estimate,conf_low,conf_high)
-
-tableOutput<-rbind(tableOutput,data.frame(
- formula="month stratified",
- VE=glue("{100 * round(VE_calc$estimate,3)}% ({100 * round(VE_calc$conf_low,3)}% - {100 * round(VE_calc$conf_high,3)}%)"),
- additional="")
-)
-
-VE_calc<-clogit(formula=case ~ vax_status + strata(siteID) + strata(month), data=vaxData)%>%
- tidy(conf.int=TRUE, exponentiate=TRUE)%>%
- mutate(estimate=1-estimate,
- conf_high=1-conf.low,
- conf_low=1-conf.high)%>%
- select(estimate,conf_low,conf_high)
-
-tableOutput<-rbind(tableOutput,data.frame(
- formula="site/month stratified",
- VE=glue("{100 * round(VE_calc$estimate,3)}% ({100 * round(VE_calc$conf_low,3)}% - {100 * round(VE_calc$conf_high,3)}%)"),
- additional="")
-)
-
-VE_calc<-clogit(formula=case ~ vax_status + strata(siteID) + strata(month) + breast, data=vaxData)%>%
- tidy(conf.int=TRUE, exponentiate=TRUE)%>%
- mutate(estimate=1-estimate,
- conf_high=1-conf.low,
- conf_low=1-conf.high)%>%
- select(estimate,conf_low,conf_high)
-
-tableOutput<-rbind(tableOutput,data.frame(
- formula="site/month stratified + BF",
- VE=glue("{100 * round(VE_calc[1,]$estimate,3)}% ({100 * round(VE_calc[1,]$conf_low,3)}% - {100 * round(VE_calc[1,]$conf_high,3)}%)"),
- additional=glue("Breastfeeding OR {1-round(VE_calc[2,]$estimate,2)} ({1-round(VE_calc[2,]$conf_high,2)} - {1-round(VE_calc[2,]$conf_low,2)})"))
-)
-
-VE_calc<-clogit(formula=case ~ vax_status + strata(siteID) + strata(month) + age_under_3_months, data=vaxData)%>%
- tidy(conf.int=TRUE, exponentiate=TRUE)%>%
- mutate(estimate=1-estimate,
- conf_high=1-conf.low,
- conf_low=1-conf.high)%>%
- select(estimate,conf_low,conf_high)
-
-tableOutput<-rbind(tableOutput,data.frame(
- formula="site/month stratified + age",
- VE=glue("{100 * round(VE_calc[1,]$estimate,3)}% ({100 * round(VE_calc[1,]$conf_low,3)}% - {100 * round(VE_calc[1,]$conf_high,3)}%)"),
- additional=glue("<3m OR {1-round(VE_calc[2,]$estimate,2)} ({1-round(VE_calc[2,]$conf_high,2)} - {1-round(VE_calc[2,]$conf_low,2)})"))
-)
-
-VE_calc<-clogit(formula=case ~ vax_status + strata(siteID) + strata(month) + age_under_3_months, data=vaxData)%>%
- tidy(conf.int=TRUE, exponentiate=TRUE)%>%
- mutate(estimate=1-estimate,
- conf_high=1-conf.low,
- conf_low=1-conf.high)%>%
- select(estimate,conf_low,conf_high)
-
-tableOutput<-rbind(tableOutput,data.frame(
- formula="site/month stratified + age",
- VE=glue("{100 * round(VE_calc[1,]$estimate,3)}% ({100 * round(VE_calc[1,]$conf_low,3)}% - {100 * round(VE_calc[1,]$conf_high,3)}%)"),
- additional=glue("<3m OR {1-round(VE_calc[2,]$estimate,2)} ({1-round(VE_calc[2,]$conf_high,2)} - {1-round(VE_calc[2,]$conf_low,2)})"))
-)
-
-
-
-VE_calc<-clogit(formula=case ~ pertussis_status, data=vaxData)%>%
- tidy(conf.int=TRUE, exponentiate=TRUE)%>%
- mutate(estimate=1-estimate,
- conf_high=1-conf.low,
- conf_low=1-conf.high)%>%
- select(estimate,conf_low,conf_high)
-
-tableOutput<-rbind(tableOutput,data.frame(
- formula="pertussis VE vs RSV",
- VE=glue("{100 * round(VE_calc$estimate,3)}% ({100 * round(VE_calc$conf_low,3)}% - {100 * round(VE_calc$conf_high,3)}%)"),
- additional="")
-)
-
-VE_calc<-clogit(formula=case ~ pertussis_status + strata(siteID), data=vaxData)%>%
- tidy(conf.int=TRUE, exponentiate=TRUE)%>%
- mutate(estimate=1-estimate,
- conf_high=1-conf.low,
- conf_low=1-conf.high)%>%
- select(estimate,conf_low,conf_high)
-
-tableOutput<-rbind(tableOutput,data.frame(
- formula="pertussis VE vs RSV (site stratified)",
- VE=glue("{100 * round(VE_calc$estimate,3)}% ({100 * round(VE_calc$conf_low,3)}% - {100 * round(VE_calc$conf_high,3)}%)"),
- additional="")
-)
-
-
-t2<-tableOutput%>%gt()
-
- return(list(table1=t1,table2=t2))
-}
-
-calculate_adjusted_VE<-function(vaxData){
-
- VE_calc<-clogit(formula=case ~ vax_status + strata(siteID) + strata(month) + age_under_3_months + male_gender + prem, data=vaxData)%>%
- tidy(conf.int=TRUE, exponentiate=TRUE)%>%
- mutate(VEestimate=1-estimate, VEconf_high=1-conf.low, VEconf_low=1-conf.high)%>%
- select(estimate,conf.low,conf.high,VEestimate,VEconf_high,VEconf_low)
-
-tableOutput<-data.frame(variable=NULL,estimate=NULL)
-tableOutput<-rbind(tableOutput,data.frame(
- variable="vaccine effectiveness",
- estimate=glue("{100 * round(VE_calc[1,]$VEestimate,3)}% ({100 * round(VE_calc[1,]$VEconf_low,3)}% - {100 * round(VE_calc[1,]$VEconf_high,3)}%)")))
-
-tableOutput<-rbind(tableOutput,data.frame(
- variable="OR age <3m",
- estimate=glue("{round(VE_calc[2,]$estimate,2)} ({round(VE_calc[2,]$conf.low,2)} - {round(VE_calc[2,]$conf.high,2)})")))
-
-tableOutput<-rbind(tableOutput,data.frame(
- variable="OR male",
- estimate=glue("{round(VE_calc[3,]$estimate,2)} ({round(VE_calc[3,]$conf.low,2)} - {round(VE_calc[3,]$conf.high,2)})")))
-
-tableOutput<-rbind(tableOutput,data.frame(
- variable="OR <37w gestation",
- estimate=glue("{round(VE_calc[4,]$estimate,2)} ({round(VE_calc[4,]$conf.low,2)} - {round(VE_calc[4,]$conf.high,2)})")))
-
-tableOutput%>%gt()%>%tab_header(title="fully adjusted model",subtitle = "case ~ vax_status + strata(siteID) + strata(month) + age_under_3_months + male_gender + prem")%>%return()
-
-}
-
-make_summary_tables<-function(vaxData){
-t1<-vaxData%>%
- drop_na(case)%>%
- mutate(case=if_else(case,"case","control"))%>%
- mutate(male_gender=if_else(male_gender,"Male","Female"))%>%
- mutate(prem=if_else(prem,"Born preterm (<37w)","Born at term"))%>%
- mutate(maternal_white=if_else(maternal_white,"White","non-white"))%>%
- mutate(breastfed=if_else(breastfed,"Yes","No"))%>%
- mutate(any_oxygen=treat_hosp_3|treat_hosp_4|treat_hosp_5|treat_hosp_6)%>%
- tbl_summary(include = c(age_mnths_calc, male_gender, prem,comorb_2:comorb_6 ,maternal_white,breastfed,imd,
- any_oxygen,treat_hosp_4,treat_hosp_5,treat_hosp_6,admit_locn_3,admit_locn_4),
- type=list(age_mnths_calc~"continuous",
- male_gender~"categorical",
- prem~"categorical",
- breastfed~"categorical",
- imd~"continuous",
- maternal_white="categorical",
- any_oxygen="dichotomous",
- treat_hosp_4="dichotomous",
- treat_hosp_5="dichotomous",
- treat_hosp_6="dichotomous",
- admit_locn_3="dichotomous",
- admit_locn_4="dichotomous"
- ),
- label=(list(age_mnths_calc~"Age in months",
- male_gender~"Male gender",
- prem~"Gestation at birth",
- breastfed ="Breastfed",
- imd="Socioeconomic quintile",
- maternal_white="Maternal ethnicity",
- comorb_2="Chronic lung disease",
- comorb_3="Congenital cardiac disease",
- comorb_4="Neuromuscular disease",
- comorb_5="Other",
- comorb_6="None",
- any_oxygen="Supplemental oxygen use",
- treat_hosp_4="High flow nasal cannulae",
- treat_hosp_5="Continuous positive airway pressure",
- treat_hosp_6="Invasive ventilation",
- admit_locn_3="HDU admission",
- admit_locn_4="PICU admission"
- )),
- by=case)%>%
- add_p()%>%
- modify_caption("**Table 1. Characteristics of the included case and control patients **")
-
-t2<-vaxData%>%
- filter(case==T)%>%
- mutate(vax_status=if_else(vax_status,"vaccinated","un-vaccinated"))%>%
- tbl_summary(include = c(LOA, treat_hosp_1:treat_hosp_6,admit_locn_3,admit_locn_4),
- label=(list(
- LOA="Length of admission",
- treat_hosp_1="NG feeds",
- treat_hosp_2="IV fluids",
- treat_hosp_3="Low flow O2",
- treat_hosp_4="HFNCO2",
- treat_hosp_5="CPAP",
- treat_hosp_6="Ventilation",
- admit_locn_3="admitted to HDU",
- admit_locn_4="admitted to PICU"
- )),
- by=vax_status)%>%
- add_p()%>%
- modify_caption("**Table 2. Comparison of clinical outcomes for RSV-positive infants with mothers vaccinated more than 14 days prior to delivery compared to infants with unvaccinated mothers **")
-
-
-
-
- return(list(table1=t1,table2=t2))
-}
-
-```
-\newpage
-# Outputs using maternal recall (whole cohort)
-```{r VE_estimates_maternal}
-#| echo: false
-#| message: false
-
-output<-data%>%
- select(case,vax_status,siteID,month,breast,age_under_3_months,age_mnths_calc,pertussis_status)%>%
- drop_na(case, vax_status)%>%
- calculate_VE()
-
-output$table1
-
-output$table2%>%tab_header(title = md("VE maternal recall"))
-
-data%>%
- select(case,vax_status,siteID,month,breast,age_under_3_months,age_mnths_calc,pertussis_status, male_gender,prem)%>%
- drop_na(case, vax_status)%>%
- calculate_adjusted_VE()
-
-tables<-data%>%make_summary_tables()
-tables$table1
-tables$table2
-
-```
-\newpage
-# VE estimates excluding those who we think had it <2 weeks before birth
-```{r VE_estimates_under2weeks1}
-#| echo: false
-#| message: false
-
-remove_under_2weeks<-data%>%filter(under2weeks %in% c(NA,FALSE))
-
-remove_under_2weeks%>% mutate(vax_status=if_else(vax_status,"vax","unvax"))%>%
- tabyl(vax_status,under2weeks)%>%rename(" "=vax_status)%>%
- gt()%>%tab_header(title = md("subjects known to be vaccinated <2 weeks before birth"))
-
-remove_under_2weeks%>%
- mutate(case=if_else(case,"case","control"))%>%
- mutate(vax_status=if_else(vax_status,"vax","unvax"))%>%
- tabyl(vax_status,case)%>%
- adorn_totals(where=c("row", "col"))%>%
- rename(" "=vax_status)%>%
- gt()%>%tab_header(title = md("2x2 data feeding into regression"))
-
-remove_under_2weeks%>%
- calculate_adjusted_VE()%>%tab_header("Fully adjusted model removing <2 weeks")
-
-tables<-remove_under_2weeks%>%make_summary_tables()
-tables$table1
-tables$table2
-
-```
-
-\newpage
-# VE estimates including only those who we know had it >2 weeks before birth
-```{r VE_estimates_under2weeks2}
-#| echo: false
-#| message: false
-
-only_over_2weeks<-data%>%filter((vax_status & !under2weeks) | !vax_status)
-
-only_over_2weeks%>% mutate(vax_status=if_else(vax_status,"vax","unvax"))%>%
- tabyl(vax_status,under2weeks)%>%rename(" "=vax_status)%>%
- gt()%>%tab_header(title = md("subjects known to be vaccinated <2 weeks before birth"))
-
-only_over_2weeks%>%
- mutate(case=if_else(case,"case","control"))%>%
- mutate(vax_status=if_else(vax_status,"vax","unvax"))%>%
- tabyl(vax_status,case)%>%
- adorn_totals(where=c("row", "col"))%>%
- rename(" "=vax_status)%>%
- gt()%>%tab_header(title = md("2x2 data feeding into regression"))
-
-output<-only_over_2weeks%>%
- calculate_VE()
-output$table1
-output$table2%>%tab_header(title = md("VE only under 2 weeks"))
-
-only_over_2weeks%>%
- calculate_adjusted_VE()%>%tab_header("Fully adjusted model for known >=2 weeks")
-
-tables<-only_over_2weeks%>%make_summary_tables()
-tables$table1
-tables$table2
-
-```
-
-\newpage
-# Fully Adjusted pertussis
-```{r VE_adj_pertussis}
-#| echo: false
-#| message: false
-
-#have to fudge month & swap vax for pertussus...
-pertussis<-data%>%
- mutate(vax_status=pertussis_status)%>%
- mutate(siteID=1)
-
-pertussis%>%
- mutate(case=if_else(case,"case","control"))%>%
- mutate(vax_status=if_else(vax_status,"vax","unvax"))%>%
- tabyl(vax_status,case)%>%
- adorn_totals(where=c("row", "col"))%>%
- rename(" "=vax_status)%>%
- gt()%>%tab_header(title = md("2x2 data feeding into regression"))
-
-pertussis%>%
- calculate_adjusted_VE()%>%tab_header("Fully adjusted model for pertussis")
-
-```
-
-\newpage
-# VE estimates using the records data we have
-```{r VE_estimates_records}
-#| echo: false
-#| message: false
-
-records_only<-data%>%mutate(vax_status = case_when(vax_yn_record=="confirmed_vax" ~ TRUE,
- vax_yn_record=="confirmed_unvax" ~ FALSE))
-
-records_only%>%
- mutate(case=if_else(case,"case","control"))%>%
- mutate(vax_status=if_else(vax_status,"vax","unvax"))%>%
- tabyl(vax_status,case)%>%
- adorn_totals(where=c("row", "col"))%>%
- rename(" "=vax_status)%>%
- gt()%>%tab_header(title = md("2x2 data feeding into regression"))
-
-records_only%>%
- calculate_adjusted_VE()%>%tab_header("Just using the maternal records we have")
-
-tables<-records_only%>%make_summary_tables()
-tables$table1
-tables$table2
-
-records_only%>%filter((vax_status & !under2weeks) | !vax_status)%>%
- mutate(case=if_else(case,"case","control"))%>%
- mutate(vax_status=if_else(vax_status,"vax","unvax"))%>%
- tabyl(vax_status,case)%>%
- adorn_totals(where=c("row", "col"))%>%
- rename(" "=vax_status)%>%
- gt()%>%tab_header(title = md("2x2 data feeding into regression if we only included <2weeks"))
-
-records_only%>%filter((vax_status & !under2weeks) | !vax_status)%>%
- calculate_adjusted_VE()%>%tab_header("Just using the maternal records we have confirmed over 2 weeks!")
-
-```
-
-```{r VE_estimates_combined}
-#| echo: false
-#| message: false
-#| include: false
-
-# VE estimates using combined maternal + any records data we have
-
-output<-data%>%
- mutate(vax_status = vax_status_combined)%>%
- select(case,vax_status,siteID,month,breast,age_under_3_months,age_mnths_calc,pertussis_status)%>%
- drop_na(case, vax_status)%>%
- calculate_VE()
-
-output$table1
-
-output$table2%>%tab_header(title = md("VE combined records data"))
-
-```
-
-\newpage
-## Screening method comparision of subgroups vs all controls (unadjusted)
-(nb I don't rank the other treatments at the moment.... so no cpap can include ventilation.)
-```{r screening_method}
-#| echo: false
-#| message: false
-
-subgroupVE<-function(subgroup,outcome){
-VE_calc<-clogit(formula=case ~ vax_status, data=subgroup)%>%
- tidy(conf.int=TRUE, exponentiate=TRUE)%>%
- mutate(estimate=1-estimate,
- conf_high=1-conf.low,
- conf_low=1-conf.high)%>%
- mutate(measure=glue("VE against {outcome}"))%>%
- select(measure,estimate,conf_low,conf_high)
-cat(glue("{VE_calc$measure} {100 * round(VE_calc$estimate,3)}% ({100 * round(VE_calc$conf_low,3)}% - {100 * round(VE_calc$conf_high,3)}%)"))
-
-subgroup%>%
- mutate(case=if_else(case,"case","control"))%>%
- mutate(vax_status=if_else(vax_status,"vax","unvax"))%>%
- tabyl(case,vax_status)%>%adorn_totals(where=c("row", "col"))%>%rename(" "=case)%>%
- mutate(percVac=vax/Total)%>%mutate(percUnvac=unvax/Total)%>%
- gt()%>%fmt_percent(starts_with("perc"),decimals = 1)%>%
- tab_header(title = glue("{outcome} vs all controls"))%>%
- return()
-}
-
-outcomes<-data%>%
- select(case,vax_status,siteID,month,breast,age_under_3_months,age_mnths_calc,pertussis_status,treat_hosp_4,treat_hosp_5,treat_hosp_6,
- admit_locn_3,admit_locn_4,under2weeks)%>%
- rename(treat_hfno=treat_hosp_4,treat_cpap=treat_hosp_5,treat_vent=treat_hosp_6,
- admit_hdu=admit_locn_3,admit_picu=admit_locn_4)%>%
- drop_na(case, vax_status)
-
-outcomes%>%filter(case & treat_hfno | !case)%>%subgroupVE("HFNCO")
-
-outcomes%>%filter(case & treat_cpap | !case)%>%subgroupVE("CPAP")
-
-outcomes%>%filter(case & treat_vent | !case)%>%subgroupVE("vent")
-
-outcomes%>%filter(case & admit_hdu | !case)%>%subgroupVE("hdu")
-
-outcomes%>%filter(case & admit_picu | !case)%>%subgroupVE("picu")
-
-```
-\newpage
-## Screening method comparision of subgroups vs all controls and only under 2 weeks (unadjusted)
-```{r screening_method_over_2_weeks}
-#| echo: false
-#| message: false
-
-outcomes<-outcomes%>%
- filter((vax_status & !under2weeks) | !vax_status)
-
-outcomes%>%filter(case & treat_hfno | !case)%>%subgroupVE("HFNCO")
-
-outcomes%>%filter(case & treat_cpap | !case)%>%subgroupVE("CPAP")
-
-outcomes%>%filter(case & treat_vent | !case)%>%subgroupVE("vent")
-
-outcomes%>%filter(case & admit_hdu | !case)%>%subgroupVE("hdu")
-
-outcomes%>%filter(case & admit_picu | !case)%>%subgroupVE("picu")
-
-```
-
-\newpage
-## Screening method comparision of subgroups vs all controls. (stratified by month and site)
-can't do more adjustments or it breaks
-```{r screening_method_stratified}
-#| echo: false
-#| message: false
-
-subgroupVE<-function(subgroup,outcome){
-VE_calc<-clogit(formula=case ~ vax_status + strata(siteID) + strata(month), data=subgroup)%>%
- tidy(conf.int=TRUE, exponentiate=TRUE)%>%
- mutate(estimate=1-estimate,
- conf_high=1-conf.low,
- conf_low=1-conf.high)%>%
- mutate(measure=glue("VE against {outcome}"))%>%
- select(measure,estimate,conf_low,conf_high)
-cat(glue("{VE_calc$measure} {100 * round(VE_calc$estimate,3)}% ({100 * round(VE_calc$conf_low,3)}% - {100 * round(VE_calc$conf_high,3)}%)"))
-
-subgroup%>%
- mutate(case=if_else(case,"case","control"))%>%
- mutate(vax_status=if_else(vax_status,"vax","unvax"))%>%
- tabyl(case,vax_status)%>%adorn_totals(where=c("row", "col"))%>%
- mutate(percVac=vax/Total)%>%mutate(percUnvac=unvax/Total)%>%
- rename(" "=case)%>%gt()%>%fmt_percent(starts_with("perc"),decimals = 1)%>%
- tab_header(title = glue("{outcome} vs all controls"))%>%
- return()
-}
-
-
-
-outcomes<-data%>%
- select(case,vax_status,siteID,month,breast,age_under_3_months,age_mnths_calc,pertussis_status,treat_hosp_4,treat_hosp_5,treat_hosp_6,
- admit_locn_3,admit_locn_4,under2weeks,male_gender,prem)%>%
- rename(treat_hfno=treat_hosp_4,treat_cpap=treat_hosp_5,treat_vent=treat_hosp_6,
- admit_hdu=admit_locn_3,admit_picu=admit_locn_4)%>%
- drop_na(case, vax_status)
-
-outcomes%>%filter(case & treat_hfno | !case)%>%subgroupVE("HFNCO")
-
-outcomes%>%filter(case & treat_cpap | !case)%>%subgroupVE("CPAP")
-
-outcomes%>%filter(case & treat_vent | !case)%>%subgroupVE("vent")
-
-outcomes%>%filter(case & admit_hdu | !case)%>%subgroupVE("hdu")
-
-outcomes%>%filter(case & admit_picu | !case)%>%subgroupVE("picu")
-
-```