############################################################################################################################################# * ### Study: G2 respiratory viruses ### ### Title: Respiratory Virus Shedding in Exhaled Breath and Efficacy of Face Masks ### ### Script title: Efficacy of surgical face masks ### ### Institute: School of Public Health, The University of Hong Kong ### ### Author: Nancy Leung ### ### Date: 2020-03-14 ### ### Note: ### ### - This script is to generate: ### ### - Figure 1, Extended Data Figure 4, Extended Data Figure 5, Extended Data Figure 9 ### ### - Table 1b, Supplementary Table 2, Supplementary Table 3 ### ############################################################################################################################################# * # Libraries ---------------------------------------------------------------------- #### library(openxlsx) library(lubridate) library(tidyverse) library(VGAM) packageVersion('VGAM') # [1] ‘1.1.1’ # Data -------------------------------------------------------------------------- #### s2 <- read.csv('G2resp_derived data_log10vl.csv', stringsAsFactors = F) dim(s2) # [1] 246 279 s2 <- s2 %>% mutate_at(vars(MMWRweek.enddate, symponset), ymd) # PCR positive variables grep('pcrgp_', names(s2), value = T) # [1] "pcrgp_hcov_nasal" "pcrgp_hcov_throat" "pcrgp_hcov_conc_nm" "pcrgp_hcov_conc_wm" "pcrgp_hcov_imp_nm" "pcrgp_hcov_imp_wm" # [7] "pcrgp_hcovnl63_nasal" "pcrgp_hcovnl63_throat" "pcrgp_hcovnl63_conc_nm" "pcrgp_hcovnl63_conc_wm" "pcrgp_hcovnl63_imp_nm" "pcrgp_hcovnl63_imp_wm" # [13] "pcrgp_hcov229e_nasal" "pcrgp_hcov229e_throat" "pcrgp_hcov229e_conc_nm" "pcrgp_hcov229e_conc_wm" "pcrgp_hcov229e_imp_nm" "pcrgp_hcov229e_imp_wm" # [19] "pcrgp_hcovoc43_nasal" "pcrgp_hcovoc43_throat" "pcrgp_hcovoc43_conc_nm" "pcrgp_hcovoc43_conc_wm" "pcrgp_hcovoc43_imp_nm" "pcrgp_hcovoc43_imp_wm" # [25] "pcrgp_hcovhku1_nasal" "pcrgp_hcovhku1_throat" "pcrgp_hcovhku1_conc_nm" "pcrgp_hcovhku1_conc_wm" "pcrgp_hcovhku1_imp_nm" "pcrgp_hcovhku1_imp_wm" # [31] "pcrgp_flu_nasal" "pcrgp_flu_throat" "pcrgp_flu_conc_nm" "pcrgp_flu_conc_wm" "pcrgp_flu_imp_nm" "pcrgp_flu_imp_wm" # [37] "pcrgp_flua_nasal" "pcrgp_flua_throat" "pcrgp_flua_conc_nm" "pcrgp_flua_conc_wm" "pcrgp_flua_imp_nm" "pcrgp_flua_imp_wm" # [43] "pcrgp_flub_nasal" "pcrgp_flub_throat" "pcrgp_flub_conc_nm" "pcrgp_flub_conc_wm" "pcrgp_flub_imp_nm" "pcrgp_flub_imp_wm" # [49] "pcrgp_enterhino_nasal" "pcrgp_enterhino_throat" "pcrgp_enterhino_conc_nm" "pcrgp_enterhino_conc_wm" "pcrgp_enterhino_imp_nm" "pcrgp_enterhino_imp_wm" # [55] "pcrgp_hmpv_nasal" "pcrgp_hmpv_throat" "pcrgp_hmpv_conc_nm" "pcrgp_hmpv_conc_wm" "pcrgp_hmpv_imp_nm" "pcrgp_hmpv_imp_wm" # [61] "pcrgp_piv_nasal" "pcrgp_piv_throat" "pcrgp_piv_conc_nm" "pcrgp_piv_conc_wm" "pcrgp_piv_imp_nm" "pcrgp_piv_imp_wm" # [67] "pcrgp_rsv_nasal" "pcrgp_rsv_throat" "pcrgp_rsv_conc_nm" "pcrgp_rsv_conc_wm" "pcrgp_rsv_imp_nm" "pcrgp_rsv_imp_wm" # [73] "pcrgp_adeno_nasal" "pcrgp_adeno_throat" "pcrgp_adeno_conc_nm" "pcrgp_adeno_conc_wm" "pcrgp_adeno_imp_nm" "pcrgp_adeno_imp_wm" # viral load variables grep('log10vlgp_', names(s2), value = T) # [1] "log10vlgp_hcov_nasal" "log10vlgp_flu_nasal" "log10vlgp_enterhino_nasal" "log10vlgp_hcov_throat" "log10vlgp_flu_throat" # [6] "log10vlgp_enterhino_throat" "log10vlgp_hcov_conc_nm" "log10vlgp_flu_conc_nm" "log10vlgp_enterhino_conc_nm" "log10vlgp_hcov_conc_wm" # [11] "log10vlgp_flu_conc_wm" "log10vlgp_enterhino_conc_wm" "log10vlgp_hcov_imp_nm" "log10vlgp_flu_imp_nm" "log10vlgp_enterhino_imp_nm" # [16] "log10vlgp_hcov_imp_wm" "log10vlgp_flu_imp_wm" "log10vlgp_enterhino_imp_wm" # Functions ----------------------------------------------------------------------- #### # function to generate table of significant testing for frequency of detection by Fisher's test #### f_fisher <- function(df, virus, sample) { # prepare dataset var.nm <- paste0('pcrgp_', virus, '_', sample, '_nm') # variable for PCR positive for without mask var.wm <- paste0('pcrgp_', virus, '_', sample, '_wm') # variable for PCR positive for with mask mask.n <- df %>% filter(Mask.wo == 1) %>% select(var.nm) %>% rename(pcrpos = var.nm) %>% mutate(mask = 0) %>% select(mask, pcrpos) mask.y <- df %>% filter(Mask.wm == 1) %>% select(var.wm) %>% rename(pcrpos = var.wm) %>% mutate(mask = 1) %>% select(mask, pcrpos) mask.ny <- rbind(mask.n, mask.y) # dataset for PCR positive for without mask and with mask # significant testing for frequency of detection if(any(mask.ny$pcrpos == 0) & any(mask.ny$pcrpos == 1)) { fisher <- table(mask.ny) %>% fisher.test() # two-sided Fisher's test pvalue.raw <- fisher$p.value } if(!(any(mask.ny$pcrpos == 0) & any(mask.ny$pcrpos == 1))) { pvalue.raw <- 1 # if no difference in outcome, p-value = 1 } # generate numbers to be put in table for frequency of detection tab <- table(mask.ny) %>% addmargins %>% data.frame %>% arrange(mask, pcrpos) num.total.collected_ctrl <- dim(mask.n)[1] num.total.collected_mask <- dim(mask.y)[1] num.total.tested_ctrl <- tab %>% filter(mask == 0 & pcrpos == 'Sum') %>% select(Freq) %>% unlist %>% unname num.total.tested_mask <- tab %>% filter(mask == 1 & pcrpos == 'Sum') %>% select(Freq) %>% unlist %>% unname num.pcrpos_ctrl <- tab %>% filter(mask == 0 & pcrpos == 1) %>% select(Freq) %>% unlist %>% unname %>% ifelse(length(.)!=0, ., 0) num.pcrpos_mask <- tab %>% filter(mask == 1 & pcrpos == 1) %>% select(Freq) %>% unlist %>% unname %>% ifelse(length(.)!=0, ., 0) percent.pcrpos_ctrl <- unname(num.pcrpos_ctrl/num.total.tested_ctrl) percent.pcrpos_mask <- unname(num.pcrpos_mask/num.total.tested_mask) # generate table for frequency of detection with p-value output <- data.frame(virus = virus, sample = sample, total.collected_ctrl = num.total.collected_ctrl, total.collected_mask = num.total.collected_mask, total.tested_ctrl = num.total.tested_ctrl, total.tested_mask = num.total.tested_mask, pcrpos_ctrl = num.pcrpos_ctrl, pcrpos_mask = num.pcrpos_mask, percent_ctrl = percent.pcrpos_ctrl, percent_mask = percent.pcrpos_mask, pvalue.raw = pvalue.raw, ctrl = paste0(num.pcrpos_ctrl, '/', num.total.tested_ctrl, ' (', specify_decimal(percent.pcrpos_ctrl, 2)*100, ')'), mask = paste0(num.pcrpos_mask, '/', num.total.tested_mask, ' (', specify_decimal(percent.pcrpos_mask, 2)*100, ')'), pvalue = ifelse(!is.na(pvalue.raw), specify_decimal(pvalue.raw, 2), NA), stringsAsFactors = F) output } # function to run Tobit regression #### f_tobit <- function(df_fortest) { # for scenarios where Tobit regression cannot be run if(length(na.omit(unique(df_fortest$logvl))) == 1 | # only one unique value for outcome (i.e. log10 virus copies per sample) length(na.omit(unique(df_fortest$predictor))) == 1 | # only one unique value for predictor (i.e. factor) length(na.omit(df_fortest$predictor)) < 3 ) # sample size <= 2 { p <- NA # p-value cannot be generated } # run Tobit regression if( !(length(na.omit(unique(df_fortest$logvl))) == 1 | length(na.omit(unique(df_fortest$predictor))) == 1 | length(na.omit(df_fortest$predictor)) < 3) ) { mod <- vglm(logvl ~ predictor, tobit(Lower = log(1.9, 10)), # equivalent to 0.30103 log10(virus copies) data = df_fortest) # Tobit model on log10 virus copies per sample summary(mod) std.error <- sqrt(diag(vcov(mod))) # standard error coeffs <- coef(mod) # coefficients confints <- cbind(ci.lower = coeffs - qnorm(0.975) * std.error, # 95% confidence intervals ci.upper = coeffs + qnorm(0.975) * std.error) pvals <- 2 * pt(abs(coef(summary(mod))[, 'z value']), df.residual(mod), lower.tail = F) # p-values # coefficient, confidence interval and p-value for predictor coeff_pred <- coeffs[3] confint_pred <- confints[3, ] pvalue_pred <- pvals[3] output <- data.frame(coeff = coeff_pred, ci.lower = confint_pred[1], ci.upper = confint_pred[2], pvalue = pvalue_pred, stringsAsFactors = F) rownames(output) <- NULL ## only needs the pvalue p <- output$pvalue p } } # function to generate plot of viral load by sample with result of significant testing for viral load by Tobit regression #### f_plot.vl.sample <- function(virus, log10vlgp_nasal, log10vlgp_throat, log10vlgp_imp_nm, log10vlgp_imp_wm, log10vlgp_conc_nm, log10vlgp_conc_wm, mtitle, xlabels, fig) { ## respiratory droplets # data vlgp2.test_imp_nm <- data.frame(mask = 0, vlgp2 = log10vlgp_imp_nm) # respiratory droplets, without mask vlgp2.test_imp_wm <- data.frame(mask = 1, vlgp2 = log10vlgp_imp_wm) # respiratory droplets, with mask fortest_imp <- rbind(vlgp2.test_imp_nm, vlgp2.test_imp_wm) %>% rename(logvl = vlgp2, predictor = mask) rownames(fortest_imp) <- NULL # significant testing by Tobit regression pvalue_imp.tobit <- f_tobit(df_fortest = fortest_imp) pvalue_imp <- ifelse(is.null(pvalue_imp.tobit), 1, pvalue_imp.tobit) ## aerosols # data vlgp2.test_conc_nm <- data.frame(mask = 0, vlgp2 = log10vlgp_conc_nm) # aerosols, without mask vlgp2.test_conc_wm <- data.frame(mask = 1, vlgp2 = log10vlgp_conc_wm) # aerosols, with mask fortest_conc <- rbind(vlgp2.test_conc_nm, vlgp2.test_conc_wm) %>% rename(logvl = vlgp2, predictor = mask) rownames(fortest_conc) <- NULL # significant testing by Tobit regression pvalue_conc.tobit <- f_tobit(df_fortest = fortest_conc) pvalue_conc <- ifelse(is.null(pvalue_conc.tobit), 1, pvalue_conc.tobit) ## plot parameters fig.xpos <- 0.24 fig.ypos <- 2.5 jitter.set <- 0.2 cex.set <- 1 cex.xlabels <- 1.3 x.axis.pos <- -0.5 xlabels1.ypos <- -1.2 xlabels.gap <- 0.5 xlabels1 <- c('Nasal', 'Throat', expression('Droplet particles' > '5 µm,'), expression('Droplet particles' > '5 µm,'), expression('Aerosol particles' <= '5 µm,'), expression('Aerosol particles' <= '5 µm,')) xlabels2 <- c('swab', 'swab', 'without mask', 'with mask', 'without mask', 'with mask') pvalue.line.y <- 8 pvalue.value.y <- pvalue.line.y + 0.5 cex.pvalue <- 1.8 plot_paras <- data.frame( sample = c('nasal', 'throat', 'imp_nm', 'imp_wm', 'conc_nm', 'conc_wm'), x.set = 1:6, col.set = c('red', 'blue', 'green4', 'mediumseagreen', 'darkorange4', 'darkorange1'), pch.set = 16, stringsAsFactors = F) x.axis.limits <- c(min(plot_paras$x.set) - 0.5, max(plot_paras$x.set) + 0.5) ## plot box plot and extract box plot output # plot bp.out <- boxplot( log10vlgp_nasal, log10vlgp_throat, log10vlgp_imp_nm, log10vlgp_imp_wm, log10vlgp_conc_nm, log10vlgp_conc_wm, at = plot_paras$x.set + 0.3, boxwex = 0.10, names = plot_paras$sample, outline = F, xlim = x.axis.limits, ylim = c(0, 10), axes = F, main = '', xlab = '', ylab = '', xpd = NA) # extract statistics from box plot output bp.out.stats <- bp.out$stats %>% data.frame colnames(bp.out.stats) <- bp.out$names rownames(bp.out.stats) <- c('min', 'lower_quartile', 'median', 'upper_quartile', 'max') bp.out.stats.t <- bp.out.stats %>% t %>% data.frame %>% mutate(sample = rownames(.), virus = virus, pvalue_imp = pvalue_imp, pvalue_conc = pvalue_conc) %>% select(virus, sample, everything()) # add x-axis and labels axis(1, at = plot_paras$x.set, labels = rep('', length(plot_paras$x.set)), lwd = 0, lwd.ticks = 1, pos = x.axis.pos) segments(x0 = x.axis.limits[1], x1 = x.axis.limits[2], y0 = x.axis.pos, xpd = NA) par(xpd = T) text(xlabels1, x = plot_paras$x.set, y = xlabels1.ypos, cex = cex.xlabels, srt = 0, adj = c(0.5, 0.5), xpd = NA) text(xlabels2, x = plot_paras$x.set, y = xlabels1.ypos - xlabels.gap, cex = cex.xlabels, srt = 0, adj = c(0.5, 0.5), xpd = NA) # add y-axis and labels yaxis.label <- c(expression(10^10), '', expression(10^8), '', expression(10^6), '', expression(10^4), '', expression(10^2), '', expression(10^0)) axis(2, at = 10:0, labels = yaxis.label, cex.axis = 1.8, las = 1, pos = 0.5) # add main title mtext(text = mtitle, side = 3, line = 1.5, cex = 1.6, font = 2) # add x-axis and y-axis titles mtext(text = "Sample type", side = 1, line = 4.8, cex = 1.3, font = 1) mtext(text = "Virus copies per sample", side = 2, line = 2, cex = 1.3, font = 1) # add figure A/B/C/D text(fig, x = 0 - fig.xpos, y = 10 + fig.ypos, cex = 3, font = 2, srt = 0, adj = c(0, 0), xpd = NA) ## add individual data points to box plot # variables for data log10vl_virus.objs <- c('log10vlgp_nasal', 'log10vlgp_throat', 'log10vlgp_imp_nm', 'log10vlgp_imp_wm', 'log10vlgp_conc_nm', 'log10vlgp_conc_wm') for (i in 1:length(log10vl_virus.objs)) { # source data log10vl_virus <- get(log10vl_virus.objs[i]) # add points set.seed(seed_num.master); points(y = log10vl_virus, x = jitter(rep(plot_paras$x.set[i], length(log10vl_virus)), amount = jitter.set), col = plot_paras$col.set[i], pch = plot_paras$pch.set[i], cex = cex.set) } ## add p-values to box plot (if applicable) # bold p-value if significant bold_imp <- ifelse(!is.na(pvalue_imp) & pvalue_imp <= 0.05, 2, 1) bold_conc <- ifelse(!is.na(pvalue_conc) & pvalue_conc <= 0.05, 2, 1) # respiratory droplets if(!is.na(pvalue_imp)) { p.value.imp.line.x1 <- plot_paras$x.set[plot_paras$sample=='imp_nm'] p.value.imp.line.x2 <- plot_paras$x.set[plot_paras$sample=='imp_wm'] segments(x0 = p.value.imp.line.x1, x1 = p.value.imp.line.x2, y0 = pvalue.line.y) text(paste0('p = ', sprintf('%.2f', specify_decimal(pvalue_imp, 2))), x = (p.value.imp.line.x1 + p.value.imp.line.x2)/2, y = pvalue.value.y, cex = cex.pvalue, font = bold_imp) } # aerosols if(!is.na(pvalue_conc)) { p.value.conc.line.x1 <- plot_paras$x.set[plot_paras$sample=='conc_nm'] p.value.conc.line.x2 <- plot_paras$x.set[plot_paras$sample=='conc_wm'] segments(x0 = p.value.conc.line.x1, x1 = p.value.conc.line.x2, y0 = pvalue.line.y) text(paste0('p = ', sprintf('%.2f', specify_decimal(pvalue_conc, 2))), x = (p.value.conc.line.x1 + p.value.conc.line.x2)/2, y = pvalue.value.y, cex = cex.pvalue, font = bold_conc) } ## export box plot output logvl_quantiles <- bp.out.stats.t logvl_quantiles } # function to generate plot of viral load by days since symptom onset #### f_plot.vl.onset <- function(virus, df, mtitle, fig) { ## generate dataset for viral load per sample per day since symptom onset as column df_s <- df[, c('pid', 'days.since.onset', paste0('log10vlgp_', virus, c('_nasal', '_throat', '_imp_nm', '_conc_nm')))] df_w.raw <- reshape(df_s, idvar = 'pid', timevar = 'days.since.onset', sep = '_d', direction = 'wide') %>% select(pid, contains('nasal'), contains('throat'), contains('imp_nm'), contains('conc_nm'), everything()) ## add place holder to dataset so that all days 0/1/2/3 has a column df_w.allnames <- paste0('log10vlgp_', virus, rep(c('_nasal', '_throat', '_imp_nm', '_conc_nm'), each = 4), '_d', rep(0:3, times = 4)) df_w.missnames <- setdiff(df_w.allnames, names(df_w.raw)) df_w_misscols <- data.frame(matrix(numeric(), dim(df_w.raw)[1], length(df_w.missnames), dimnames = list(c(), df_w.missnames)), stringsAsFactors = F) df_w <- cbind(df_w.raw, df_w_misscols) %>% select(pid, contains('_d0'), contains('_d1'), contains('_d2'), contains('_d3'), everything())%>% select(pid, contains('nasal'), contains('throat'), contains('imp_nm'), contains('conc_nm'), everything()) ## plot parameters fig.xpos <- 1.4 fig.ypos <- 3.7 jitter.set <- 0.2 cex.set <- 1 cex.xlabels <- 1.3 x.axis.pos <- -0.5 xlabels1.ypos <- -1.45 x.set.gap.within_gp <- 0.9 x.set.gap.between_gp <- 1.5 x.set.1st <- 1 + x.set.gap.within_gp*0:3 x.set.2nd <- max(x.set.1st) + x.set.gap.between_gp + x.set.gap.within_gp*0:3 x.set.3rd <- max(x.set.2nd) + x.set.gap.between_gp + x.set.gap.within_gp*0:3 x.set.4th <- max(x.set.3rd) + x.set.gap.between_gp + x.set.gap.within_gp*0:3 x.set.all <- c(x.set.1st, x.set.2nd, x.set.3rd, x.set.4th) x.set.all plot_paras <- data.frame( sample = rep(c('nasal', 'throat', 'imp_nm', 'conc_nm'), each = 4), days.since.onset = rep(0:3, times = 4), var = setdiff(names(df_w), 'pid'), x.set = x.set.all, xlabels = rep(0:3, times = 4), col.set = rep(c('red', 'blue', 'green4', 'darkorange4'), each = 4), pch.set = 16, stringsAsFactors = F) samlabels1.xpos <- c(mean(plot_paras$x.set[2:3]), mean(plot_paras$x.set[6:7]), mean(plot_paras$x.set[10:11]), mean(plot_paras$x.set[14:15])) samlabels.line.x1 <- plot_paras$x.set[plot_paras$days.since.onset==0] samlabels.line.x2 <- plot_paras$x.set[plot_paras$days.since.onset==3] samlabels1.ypos <- 11.7 samlabels.gap1 <- 0.5 samlabels.gap2 <- 0.5 cex.samlabels <- 1.5 samlabels1 <- c('Nasal', 'Throat', expression('Droplet particles' > '5 µm,'), expression('Aerosol particles' <= '5 µm,')) samlabels2 <- c('swab', 'swab', 'without mask', 'without mask') x.axis.limits <- c(min(plot_paras$x.set) - 0.5, max(plot_paras$x.set) + 0.5) ## plot box plot and extract box plot output # plot bp.out <- boxplot( df_w[plot_paras$var], at = plot_paras$x.set + 0.4, boxwex = 0.13, names = paste0(plot_paras$sample, '_d', plot_paras$days.since.onset), outline = F, xlim = x.axis.limits, ylim = c(0, 10), axes = F, main = '', xlab = '', ylab = '', xpd = NA) # extract statistics from box plot output bp.out.stats <- bp.out$stats %>% data.frame colnames(bp.out.stats) <- bp.out$names rownames(bp.out.stats) <- c('min', 'lower_quartile', 'median', 'upper_quartile', 'max') bp.out.stats.t <- bp.out.stats %>% t %>% data.frame %>% mutate(var = rownames(.), virus = virus, sample = gsub('^(.*)(_d)(.*)$', '\\1', var), days.since.onset = gsub('^(.*)(_d)(.*)$', '\\3', var) ) %>% select(virus, sample, days.since.onset, var, everything()) # add x-axis and labels axis(1, at = plot_paras$x.set, labels = rep('', dim(plot_paras)[1]), lwd = 0, lwd.ticks = 1, pos = x.axis.pos) segments(x0 = x.axis.limits[1], x1 = x.axis.limits[2], y0 = x.axis.pos, xpd = NA) par(xpd = T) text(plot_paras$xlabels, x = plot_paras$x.set, y = xlabels1.ypos, cex = cex.xlabels, srt = 0, adj = c(0.5, 0.5), xpd = NA) # add y-axis and labels yaxis.label <- c(expression(10^10), '', expression(10^8), '', expression(10^6), '', expression(10^4), '', expression(10^2), '', expression(10^0)) axis(2, at = 10:0, labels = yaxis.label, cex.axis = 1.8, las = 1, pos = 0.5) # add main title mtext(text = mtitle, side = 3, line = 4.4, cex = 1.6, font = 2) # add sample labels par(xpd = T) text(samlabels1, x = samlabels1.xpos, y = samlabels1.ypos, cex = cex.samlabels, srt = 0, adj = c(0.5, 0.5), xpd = NA) text(samlabels2, x = samlabels1.xpos, y = samlabels1.ypos - samlabels.gap1, cex = cex.samlabels, srt = 0, adj = c(0.5, 0.5), xpd = NA) for (i in 1:length(samlabels.line.x1)) { segments(x0 = samlabels.line.x1[i] - (x.set.gap.within_gp/2), x1 = samlabels.line.x2[i] + (x.set.gap.within_gp/2), y0 = samlabels1.ypos - samlabels.gap1 - samlabels.gap2) } # add x-axis and y-axis titles mtext(text = 'Days since symptom onset', side = 1, line = 4.8, cex = 1.3, font = 1) mtext(text = 'Virus copies per sample', side = 2, line = 2, cex = 1.3, font = 1) # add figure A/B/C/D text(fig, x = 0 - fig.xpos, y = 10 + fig.ypos, cex = 3, font = 2, srt = 0, adj = c(0, 0), xpd = NA) ## add individual data points to box plot for (i in 1:dim(plot_paras)[1]) { # source data vl2_virus <- df_w[plot_paras$var[i]] %>% unlist # add points set.seed(seed_num.master); points(y = vl2_virus, x = jitter(rep(plot_paras$x.set[i], length(vl2_virus)), amount = jitter.set), col = plot_paras$col.set[i], pch = plot_paras$pch.set[1], cex = cex.set) } ## export box plot output logvl_quantiles_all <- bp.out.stats.t logvl_quantiles_all } # function to bind results of frequency of detection and viral load together into one table #### f_bind_freq_vl <- function(tab_freq, tab_logvl) { # rename columns for frequency of detection tab_freq2 <- tab_freq %>% mutate(type = 'freq') %>% rename_all(~gsub('\\.', '_', .)) %>% select(virus, type, ctrl_imp, mask_imp, pvalue_imp, ctrl_conc, mask_conc, pvalue_conc) # rename columns for viral load tab_logvl2 <- tab_logvl %>% mutate(type = 'logvl') %>% rename(nasal = median_iqr.nasal, throat = median_iqr.throat, ctrl_imp = median_iqr.imp_nm, mask_imp = median_iqr.imp_wm, ctrl_conc = median_iqr.conc_nm, mask_conc = median_iqr.conc_wm) %>% select(virus, type, ctrl_imp, mask_imp, pvalue_imp, ctrl_conc, mask_conc, pvalue_conc) tab <- rbind(tab_freq2, tab_logvl2) tab } # Generate plot: Viral load by sample for coronavirus, influenza virus and rhinovirus ------------------------------------- #### ## n # number of all infected individuals s2 %>% filter(pcrgpany_hcov == 1) %>% dim %>% .[1] # [1] 17 s2 %>% filter(pcrgpany_flu == 1) %>% dim %>% .[1] # [1] 43 s2 %>% filter(pcrgpany_enterhino == 1) %>% dim %>% .[1] # [1] 54 # number of infected individuals who provided exhaled breath without wearing mask s2 %>% filter(pcrgpany_hcov == 1) %>% filter(Mask.wo == 1) %>% dim() %>% .[1] # [1] 10 s2 %>% filter(pcrgpany_flu == 1) %>% filter(Mask.wo == 1) %>% dim() %>% .[1] # [1] 23 s2 %>% filter(pcrgpany_enterhino == 1) %>% filter(Mask.wo == 1) %>% dim() %>% .[1] # [1] 36 # number of infected individuals who provided exhaled breath wearing a mask s2 %>% filter(pcrgpany_hcov == 1) %>% filter(Mask.wm == 1) %>% dim() %>% .[1] # [1] 11 s2 %>% filter(pcrgpany_flu == 1) %>% filter(Mask.wm == 1) %>% dim() %>% .[1] # [1] 28 s2 %>% filter(pcrgpany_enterhino == 1) %>% filter(Mask.wm == 1) %>% dim() %>% .[1] # [1] 32 # number of infected individuals who provided both exhaled breath without wearing and while wearing a mask s2 %>% filter(pcrgpany_hcov == 1) %>% filter(Mask.wo == 1 & Mask.wm == 1) %>% dim() %>% .[1] # [1] 4 s2 %>% filter(pcrgpany_flu == 1) %>% filter(Mask.wo == 1 & Mask.wm == 1) %>% dim() %>% .[1] # [1] 8 s2 %>% filter(pcrgpany_enterhino == 1) %>% filter(Mask.wo == 1 & Mask.wm == 1) %>% dim() %>% .[1] # [1] 14 ## plot setEPS() postscript('Figure 1.eps', width = 12, height = 15) par(mfrow = c(3, 1)) par(oma = c(2, 7, 1, 0), mar = c(8, 2, 7, 2)) # Coronavirus logvl.sample_quantiles_hcov <- f_plot.vl.sample(virus = 'hcov', log10vlgp_nasal = s2 %>% filter(pcrgpany_hcov == 1) %>% select(log10vlgp_hcov_nasal) %>% unlist, log10vlgp_throat = s2 %>% filter(pcrgpany_hcov == 1) %>% select(log10vlgp_hcov_throat) %>% unlist, log10vlgp_imp_nm = s2 %>% filter(pcrgpany_hcov == 1) %>% select(log10vlgp_hcov_imp_nm) %>% unlist, log10vlgp_imp_wm = s2 %>% filter(pcrgpany_hcov == 1) %>% select(log10vlgp_hcov_imp_wm) %>% unlist, log10vlgp_conc_nm = s2 %>% filter(pcrgpany_hcov == 1) %>% select(log10vlgp_hcov_conc_nm) %>% unlist, log10vlgp_conc_wm = s2 %>% filter(pcrgpany_hcov == 1) %>% select(log10vlgp_hcov_conc_wm) %>% unlist, mtitle = 'Coronavirus', fig = 'a' ) # Influenza virus logvl.sample_quantiles_flu <- f_plot.vl.sample(virus = 'flu', log10vlgp_nasal = s2 %>% filter(pcrgpany_flu == 1) %>% select(log10vlgp_flu_nasal) %>% unlist, log10vlgp_throat = s2 %>% filter(pcrgpany_flu == 1) %>% select(log10vlgp_flu_throat) %>% unlist, log10vlgp_imp_nm = s2 %>% filter(pcrgpany_flu == 1) %>% select(log10vlgp_flu_imp_nm) %>% unlist, log10vlgp_imp_wm = s2 %>% filter(pcrgpany_flu == 1) %>% select(log10vlgp_flu_imp_wm) %>% unlist, log10vlgp_conc_nm = s2 %>% filter(pcrgpany_flu == 1) %>% select(log10vlgp_flu_conc_nm) %>% unlist, log10vlgp_conc_wm = s2 %>% filter(pcrgpany_flu == 1) %>% select(log10vlgp_flu_conc_wm) %>% unlist, mtitle = 'Influenza virus', fig = 'b' ) # Rhinovirus logvl.sample_quantiles_enterhino <- f_plot.vl.sample(virus = 'enterhino', log10vlgp_nasal = s2 %>% filter(pcrgpany_enterhino == 1) %>% select(log10vlgp_enterhino_nasal) %>% unlist, log10vlgp_throat = s2 %>% filter(pcrgpany_enterhino == 1) %>% select(log10vlgp_enterhino_throat) %>% unlist, log10vlgp_imp_nm = s2 %>% filter(pcrgpany_enterhino == 1) %>% select(log10vlgp_enterhino_imp_nm) %>% unlist, log10vlgp_imp_wm = s2 %>% filter(pcrgpany_enterhino == 1) %>% select(log10vlgp_enterhino_imp_wm) %>% unlist, log10vlgp_conc_nm = s2 %>% filter(pcrgpany_enterhino == 1) %>% select(log10vlgp_enterhino_conc_nm) %>% unlist, log10vlgp_conc_wm = s2 %>% filter(pcrgpany_enterhino == 1) %>% select(log10vlgp_enterhino_conc_wm) %>% unlist, mtitle = 'Rhinovirus', fig = 'c' ) dev.off() # Generate plot: Viral load by sample for coronavirus NL63, OC43, HKU1, influenza A and B virus ------------------------------------- #### # Coronavirus NL63, OC43 and HKU1 #### ## n ## number of all infected individuals s2 %>% filter(pcrgpany_hcovnl63 == 1) %>% dim() %>% .[1] # [1] 8 s2 %>% filter(pcrgpany_hcovoc43 == 1) %>% dim() %>% .[1] # [1] 5 s2 %>% filter(pcrgpany_hcovhku1 == 1) %>% dim() %>% .[1] # [1] 4 ## number of infected individuals who provided exhaled breath without wearing mask s2 %>% filter(pcrgpany_hcovnl63 == 1) %>% filter(Mask.wo == 1) %>% dim() %>% .[1] # [1] 3 s2 %>% filter(pcrgpany_hcovoc43 == 1) %>% filter(Mask.wo == 1) %>% dim() %>% .[1] # [1] 3 s2 %>% filter(pcrgpany_hcovhku1 == 1) %>% filter(Mask.wo == 1) %>% dim() %>% .[1] # [1] 4 ## number of infected individuals who provided exhaled breath wearing a mask s2 %>% filter(pcrgpany_hcovnl63 == 1) %>% filter(Mask.wm == 1) %>% dim() %>% .[1] # [1] 5 s2 %>% filter(pcrgpany_hcovoc43 == 1) %>% filter(Mask.wm == 1) %>% dim() %>% .[1] # [1] 4 s2 %>% filter(pcrgpany_hcovhku1 == 1) %>% filter(Mask.wm == 1) %>% dim() %>% .[1] # [1] 2 ## number of infected individuals who provided both exhaled breath without wearing and while wearing a mask s2 %>% filter(pcrgpany_hcovnl63 == 1) %>% filter(Mask.wo == 1 & Mask.wm == 1) %>% dim() %>% .[1] # [1] 0 s2 %>% filter(pcrgpany_hcovoc43 == 1) %>% filter(Mask.wo == 1 & Mask.wm == 1) %>% dim() %>% .[1] # [1] 2 s2 %>% filter(pcrgpany_hcovhku1 == 1) %>% filter(Mask.wo == 1 & Mask.wm == 1) %>% dim() %>% .[1] # [1] 2 ## plot setEPS() postscript('ED Figure 4.eps', width = 12, height = 15) par(mfrow = c(3, 1)) par(oma = c(2, 7, 1, 0), mar = c(8, 2, 7, 2)) # coronavirus NL63 logvl_quantiles_hcovnl63 <- f_plot.vl.sample( virus = 'hcovnl63', log10vlgp_nasal = s2 %>% filter(pcrgpany_hcovnl63 == 1) %>% select(log10vl_hcovnl63_nasal) %>% unlist, log10vlgp_throat = s2 %>% filter(pcrgpany_hcovnl63 == 1) %>% select(log10vl_hcovnl63_throat) %>% unlist, log10vlgp_imp_nm = s2 %>% filter(pcrgpany_hcovnl63 == 1) %>% select(log10vl_hcovnl63_imp_nm) %>% unlist, log10vlgp_imp_wm = s2 %>% filter(pcrgpany_hcovnl63 == 1) %>% select(log10vl_hcovnl63_imp_wm) %>% unlist, log10vlgp_conc_nm = s2 %>% filter(pcrgpany_hcovnl63 == 1) %>% select(log10vl_hcovnl63_conc_nm) %>% unlist, log10vlgp_conc_wm = s2 %>% filter(pcrgpany_hcovnl63 == 1) %>% select(log10vl_hcovnl63_conc_wm) %>% unlist, mtitle = 'Coronavirus NL63', fig = 'a' ) # coronavirus OC43 logvl_quantiles_hcovoc43 <- f_plot.vl.sample( virus = 'hcovoc43', log10vlgp_nasal = s2 %>% filter(pcrgpany_hcovoc43 == 1) %>% select(log10vl_hcovoc43_nasal) %>% unlist, log10vlgp_throat = s2 %>% filter(pcrgpany_hcovoc43 == 1) %>% select(log10vl_hcovoc43_throat) %>% unlist, log10vlgp_imp_nm = s2 %>% filter(pcrgpany_hcovoc43 == 1) %>% select(log10vl_hcovoc43_imp_nm) %>% unlist, log10vlgp_imp_wm = s2 %>% filter(pcrgpany_hcovoc43 == 1) %>% select(log10vl_hcovoc43_imp_wm) %>% unlist, log10vlgp_conc_nm = s2 %>% filter(pcrgpany_hcovoc43 == 1) %>% select(log10vl_hcovoc43_conc_nm) %>% unlist, log10vlgp_conc_wm = s2 %>% filter(pcrgpany_hcovoc43 == 1) %>% select(log10vl_hcovoc43_conc_wm) %>% unlist, mtitle = 'Coronavirus OC43', fig = 'b' ) # coronavirus HKU1 logvl_quantiles_hcovhku1 <- f_plot.vl.sample( virus = 'hcovhku1', log10vlgp_nasal = s2 %>% filter(pcrgpany_hcovhku1 == 1) %>% select(log10vl_hcovhku1_nasal) %>% unlist, log10vlgp_throat = s2 %>% filter(pcrgpany_hcovhku1 == 1) %>% select(log10vl_hcovhku1_throat) %>% unlist, log10vlgp_imp_nm = s2 %>% filter(pcrgpany_hcovhku1 == 1) %>% select(log10vl_hcovhku1_imp_nm) %>% unlist, log10vlgp_imp_wm = s2 %>% filter(pcrgpany_hcovhku1 == 1) %>% select(log10vl_hcovhku1_imp_wm) %>% unlist, log10vlgp_conc_nm = s2 %>% filter(pcrgpany_hcovhku1 == 1) %>% select(log10vl_hcovhku1_conc_nm) %>% unlist, log10vlgp_conc_wm = s2 %>% filter(pcrgpany_hcovhku1 == 1) %>% select(log10vl_hcovhku1_conc_wm) %>% unlist, mtitle = 'Coronavirus HKU1', fig = 'c' ) dev.off() ## Influenza A and B virus #### ## n ## number of all infected individuals s2 %>% filter(pcrgpany_flua == 1) %>% dim() %>% .[1] # [1] 31 s2 %>% filter(pcrgpany_flub == 1) %>% dim() %>% .[1] # [1] 14 ## number of infected individuals who provided exhaled breath without wearing mask s2 %>% filter(pcrgpany_flua == 1) %>% filter(Mask.wo == 1) %>% dim() %>% .[1] # [1] 19 s2 %>% filter(pcrgpany_flub == 1) %>% filter(Mask.wo == 1) %>% dim() %>% .[1] # [1] 6 ## number of infected individuals who provided exhaled breath wearing a mask s2 %>% filter(pcrgpany_flua == 1) %>% filter(Mask.wm == 1) %>% dim() %>% .[1] # [1] 19 s2 %>% filter(pcrgpany_flub == 1) %>% filter(Mask.wm == 1) %>% dim() %>% .[1] # [1] 10 ## number of infected individuals who provided both exhaled breath without wearing and while wearing a mask s2 %>% filter(pcrgpany_flua == 1) %>% filter(Mask.wo == 1 & Mask.wm == 1) %>% dim() %>% .[1] # [1] 7 s2 %>% filter(pcrgpany_flub == 1) %>% filter(Mask.wo == 1 & Mask.wm == 1) %>% dim() %>% .[1] # [1] 2 ## plot setEPS() postscript('ED Figure 5.eps', width = 12, height = 15) par(mfrow = c(3, 1)) par(oma = c(2, 7, 1, 0), mar = c(8, 2, 7, 2)) # Influenza A virus logvl_quantiles_flua <- f_plot.vl.sample( virus = 'flua', log10vlgp_nasal = s2 %>% filter(pcrgpany_flua == 1) %>% select(log10vl_flua_nasal) %>% unlist, log10vlgp_throat = s2 %>% filter(pcrgpany_flua == 1) %>% select(log10vl_flua_throat) %>% unlist, log10vlgp_imp_nm = s2 %>% filter(pcrgpany_flua == 1) %>% select(log10vl_flua_imp_nm) %>% unlist, log10vlgp_imp_wm = s2 %>% filter(pcrgpany_flua == 1) %>% select(log10vl_flua_imp_wm) %>% unlist, log10vlgp_conc_nm = s2 %>% filter(pcrgpany_flua == 1) %>% select(log10vl_flua_conc_nm) %>% unlist, log10vlgp_conc_wm = s2 %>% filter(pcrgpany_flua == 1) %>% select(log10vl_flua_conc_wm) %>% unlist, mtitle = 'Influenza A virus', fig = 'a' ) # Influenza B virus logvl_quantiles_flub <- f_plot.vl.sample( virus = 'flub', log10vlgp_nasal = s2 %>% filter(pcrgpany_flub == 1) %>% select(log10vl_flub_nasal) %>% unlist, log10vlgp_throat = s2 %>% filter(pcrgpany_flub == 1) %>% select(log10vl_flub_throat) %>% unlist, log10vlgp_imp_nm = s2 %>% filter(pcrgpany_flub == 1) %>% select(log10vl_flub_imp_nm) %>% unlist, log10vlgp_imp_wm = s2 %>% filter(pcrgpany_flub == 1) %>% select(log10vl_flub_imp_wm) %>% unlist, log10vlgp_conc_nm = s2 %>% filter(pcrgpany_flub == 1) %>% select(log10vl_flub_conc_nm) %>% unlist, log10vlgp_conc_wm = s2 %>% filter(pcrgpany_flub == 1) %>% select(log10vl_flub_conc_wm) %>% unlist, mtitle = 'Influenza B virus', fig = 'b' ) dev.off() # Generate table: Efficacy of face masks for coronavirus, influenza virus and rhinovirus ------------------------------------- #### # - Note: Significant testing for frequency of detection and viral load for coronavirus, influenza virus and rhinovirus # frequency of detection #### e_virusgp <- rbind( f_fisher(df = s2 %>% filter(pcrgpany_hcov == 1), virus = 'hcov', sample = 'imp'), f_fisher(df = s2 %>% filter(pcrgpany_hcov == 1), virus = 'hcov', sample = 'conc'), f_fisher(df = s2 %>% filter(pcrgpany_flu == 1), virus = 'flu', sample = 'imp'), f_fisher(df = s2 %>% filter(pcrgpany_flu == 1), virus = 'flu', sample = 'conc'), f_fisher(df = s2 %>% filter(pcrgpany_enterhino == 1), virus = 'enterhino', sample = 'imp'), f_fisher(df = s2 %>% filter(pcrgpany_enterhino == 1), virus = 'enterhino', sample = 'conc') ) freq_efficacy_virusgp.print <- e_virusgp %>% select(virus, sample, ctrl, mask, pvalue) %>% reshape(., idvar = 'virus', timevar = 'sample', direction = 'wide') rownames(freq_efficacy_virusgp.print) <- NULL # viral load #### logvl_quantiles_virusgp <- rbind(logvl.sample_quantiles_hcov, logvl.sample_quantiles_flu, logvl.sample_quantiles_enterhino) logvl_quantiles_virusgp.print <- logvl_quantiles_virusgp %>% mutate_at(vars(min:max), ~specify_decimal(., 1)) %>% mutate_at(vars(pvalue_imp, pvalue_conc), ~specify_decimal(., 2)) %>% mutate(median_iqr = paste0(median, ' (', lower_quartile, ', ', upper_quartile, ')')) %>% select(sample, virus, median_iqr, pvalue_imp, pvalue_conc) %>% reshape(., idvar = 'virus', timevar = 'sample', direction = 'wide') %>% rename(pvalue_imp = pvalue_imp.imp_wm, pvalue_conc = pvalue_conc.conc_wm) %>% select( -pvalue_imp.nasal, -pvalue_imp.throat, -pvalue_imp.imp_nm, -pvalue_imp.conc_nm, -pvalue_imp.conc_wm, -pvalue_conc.nasal, -pvalue_conc.throat, -pvalue_conc.imp_wm, -pvalue_conc.imp_nm, -pvalue_conc.conc_nm) rownames(logvl_quantiles_virusgp.print) <- NULL # consolidated table: frequency of detection and viral load #### efficacy_virusgp <- f_bind_freq_vl(tab_freq = freq_efficacy_virusgp.print, tab_logvl = logvl_quantiles_virusgp.print) # export write.csv(efficacy_virusgp, 'Table 1b.csv', row.names = F) # Efficacy of face masks for coronavirus, influenza virus and rhinovirus # Generate table: Efficacy of face masks for coronavirus NL63, OC43, HKU1, influenza A and B virus ------------------------------------- #### # - Note: Significant testing for frequency of detection and viral load for coronavirus NL63, OC43, HKU1, influenza A and B virus # frequency of detection #### e_ind <- rbind( f_fisher(df = s2 %>% filter(pcrgpany_hcovnl63 == 1), virus = 'hcovnl63', sample = 'imp'), f_fisher(df = s2 %>% filter(pcrgpany_hcovnl63 == 1), virus = 'hcovnl63', sample = 'conc'), f_fisher(df = s2 %>% filter(pcrgpany_hcovoc43 == 1), virus = 'hcovoc43', sample = 'imp'), f_fisher(df = s2 %>% filter(pcrgpany_hcovoc43 == 1), virus = 'hcovoc43', sample = 'conc'), f_fisher(df = s2 %>% filter(pcrgpany_hcovhku1 == 1), virus = 'hcovhku1', sample = 'imp'), f_fisher(df = s2 %>% filter(pcrgpany_hcovhku1 == 1), virus = 'hcovhku1', sample = 'conc'), f_fisher(df = s2 %>% filter(pcrgpany_flua == 1), virus = 'flua', sample = 'imp'), f_fisher(df = s2 %>% filter(pcrgpany_flua == 1), virus = 'flua', sample = 'conc'), f_fisher(df = s2 %>% filter(pcrgpany_flub == 1), virus = 'flub', sample = 'imp'), f_fisher(df = s2 %>% filter(pcrgpany_flub == 1), virus = 'flub', sample = 'conc')) freq_efficacy_ind.print <- e_ind %>% select(virus, sample, ctrl, mask, pvalue) %>% reshape(., idvar = 'virus', timevar = 'sample', direction = 'wide') rownames(freq_efficacy_ind.print) <- NULL # viral load #### # coronavirus NL63, OC43, HKU1 logvl_quantiles_hcovs <- rbind(logvl_quantiles_hcovnl63, logvl_quantiles_hcovoc43, logvl_quantiles_hcovhku1) logvl_quantiles_hcovs.print <- logvl_quantiles_hcovs %>% mutate_at(vars(min:max), ~specify_decimal(., 1)) %>% mutate_at(vars(pvalue_imp, pvalue_conc), ~specify_decimal(., 2)) %>% mutate(median_iqr = paste0(median, ' (', lower_quartile, ', ', upper_quartile, ')')) %>% select(sample, virus, median_iqr, pvalue_imp, pvalue_conc) %>% reshape(., idvar = 'virus', timevar = 'sample', direction = 'wide') %>% rename(pvalue_imp = pvalue_imp.imp_wm, pvalue_conc = pvalue_conc.conc_wm) %>% select( -pvalue_imp.nasal, -pvalue_imp.throat, -pvalue_imp.imp_nm, -pvalue_imp.conc_nm, -pvalue_imp.conc_wm, -pvalue_conc.nasal, -pvalue_conc.throat, -pvalue_conc.imp_wm, -pvalue_conc.imp_nm, -pvalue_conc.conc_nm) row.names(logvl_quantiles_hcovs.print) <- NULL # influenza A and B virus logvl_quantiles_fluab <- rbind(logvl_quantiles_flua, logvl_quantiles_flub) logvl_quantiles_fluab.print <- logvl_quantiles_fluab %>% mutate_at(vars(min:max), ~specify_decimal(., 1)) %>% mutate_at(vars(pvalue_imp, pvalue_conc), ~specify_decimal(., 2)) %>% mutate(median_iqr = paste0(median, ' (', lower_quartile, ', ', upper_quartile, ')')) %>% select(sample, virus, median_iqr, pvalue_imp, pvalue_conc) %>% reshape(., idvar = 'virus', timevar = 'sample', direction = 'wide') %>% rename(pvalue_imp = pvalue_imp.imp_wm, pvalue_conc = pvalue_conc.conc_wm) %>% select( -pvalue_imp.nasal, -pvalue_imp.throat, -pvalue_imp.imp_nm, -pvalue_imp.conc_nm, -pvalue_imp.conc_wm, -pvalue_conc.nasal, -pvalue_conc.throat, -pvalue_conc.imp_wm, -pvalue_conc.imp_nm, -pvalue_conc.conc_nm) rownames(logvl_quantiles_fluab.print) <- NULL # bind viral load for coronavirus NL63, OC43, HKU1 and iinfluenza A and B virus together into one table logvl_quantiles_ind.print <- rbind(logvl_quantiles_hcovs.print, logvl_quantiles_fluab.print) # consolidated table: frequency of detection and viral load #### efficacy_ind <- f_bind_freq_vl(tab_freq = freq_efficacy_ind.print, tab_logvl = logvl_quantiles_ind.print) # export write.csv(efficacy_ind, 'Supplementary Table 2.csv', row.names = F) # Efficacy of face masks for coronavirus NL63, OC43, HKU1, influenza A and B virus # Generate plot: Viral load by days since symptom onset for coronavirus, influenza virus and rhinovirus ------------------------------------- #### setEPS() postscript('ED Figure 9.eps', width = 12, height = 15) par(mfrow = c(3, 1)) par(oma = c(0, 7, 3, 0), mar = c(8, 2, 7, 2)) # Coronavirus logvl.onset_quantiles_hcov <- f_plot.vl.onset(virus = 'hcov', df = s2 %>% filter(pcrgpany_hcov == 1), mtitle = 'Coronavirus', fig = 'a') # Influenza virus logvl.onset_quantiles_flu <- f_plot.vl.onset(virus = 'flu', df = s2 %>% filter(pcrgpany_flu == 1), mtitle = 'Influenza virus', fig = 'b') # Rhinovirus logvl.onset_quantiles_enterhino <- f_plot.vl.onset(virus = 'enterhino', df = s2 %>% filter(pcrgpany_enterhino == 1), mtitle = 'Rhinovirus', fig = 'c') dev.off() # Generate table: Viral load by days since symptom onset for coronavirus, influenza virus and rhinovirus ------------------------------------- #### logvl_quantiles_onset <- rbind(logvl.onset_quantiles_hcov, logvl.onset_quantiles_flu, logvl.onset_quantiles_enterhino) logvl_quantiles_onset.print <- logvl_quantiles_onset %>% mutate_at(vars(min:max), ~specify_decimal(., 1)) %>% mutate(median_iqr = paste0(median, ' (', lower_quartile, ', ', upper_quartile, ')')) %>% select(virus, sample, days.since.onset, median_iqr) %>% reshape(., idvar = c('virus', 'days.since.onset'), timevar = 'sample', direction = 'wide') %>% select(virus, days.since.onset, median_iqr.nasal, median_iqr.throat, median_iqr.imp_nm, median_iqr.conc_nm) rownames(logvl_quantiles_onset.print) <- NULL # export write.csv(logvl_quantiles_onset.print, 'Supplementary Table 3.csv', row.names = F) # Viral load by days since symptom onset for coronavirus, influenza virus and rhinovirus # Final clean up ---------------------------------------------------------- #### ### Remove objects not in master script rm(list=setdiff(ls(), master.objects)) ####### End of script #######