Infant ages

exposure_data %>%
  select(FamilyID, enrollment_age) %>%
  distinct() %>%
  summarize(
        N = n(),
        enroll_median_age_days = median(enrollment_age),
        IQR = paste(quantile(enrollment_age, c(0.25, 0.75)), collapse = ", "),
        range_days = paste(range(enrollment_age), collapse = ", ")
        ) %>%
  kable() %>% kable_styling(full_width = F)

Mom HIV

exposure_data %>% select(FamilyID, momhiv) %>% 
  distinct() %>%
  group_by(momhiv) %>%
  summarize(N = n()) %>% 
  kable() %>% kable_styling(full_width = F)

Overall

exposure_data_long %>% 
  group_by(virus, idpar) %>% 
  summarize(
    total = n(),
    total_observed = total - sum(interpolated),
    total_interpolate = stat_paste(sum(interpolated),  100*mean(interpolated), digits = 1)
  ) %>%
  kable() %>%
  kable_styling(full_width = F)

exposure_data_summary = exposure_data_long %>%
  mutate(
    pos_count = count > 0
    ) %>%
  subset(!interpolated) %>%
  group_by(virus, FamilyID, obs_infected, idpar) %>%
  mutate(
    total_pos = sum(pos_count),
    pct_pos = 100 * mean(pos_count)
    ) %>%
  group_by(pct_pos, total_pos, add = T) %>%
  summarise_at(vars(count), list(~n(), mean = mean, median = median, maximum = max)) %>%
  rename(N = n)

plot_labels = exposure_data_summary %>%
  gather(stat, estimate, mean, maximum, pct_pos) %>%
  group_by(stat) %>%
  summarize(min_lim = min(estimate), max_lim = ceiling(max(estimate))) %>%
  left_join(tibble(stat = c("pct_pos", "mean", "maximum"), 
                   out_lab = c("Percent~Positive", "Mean~Log[10]~VL", "Max.~Log[10]~VL"))) %>%
  mutate(stat = factor(stat, levels = c("pct_pos", "mean", "maximum"))) %>%
  arrange(stat) %>%
  ungroup() %>%
  mutate(letter_code = 1:3)

Joining, by = "stat"

pl_exposure = map(plot_labels$stat %>% levels(), function(s){

    tmp_theme = theme(
      legend.position = "none",
      axis.title = element_text(size = 10),
      axis.text = element_text(size = 9)
    )
    
    pl_lab = subset(plot_labels, stat == s)
    out_label = pl_lab$out_lab[1]
    lower_limit = pl_lab$min_lim
    upper_limit = pl_lab$max_lim
    
    pl1 = exposure_data_summary %>%
      gather(stat, estimate, mean, maximum, pct_pos) %>%
      select(-median, -total_pos, -N) %>%
      spread(idpar, estimate) %>%
      filter(stat == s) %>%
      ggplot(aes(x = S, y = M, colour = factor(obs_infected))) +
      geom_point() +
      geom_abline() +
      geom_point() +
      scale_colour_manual(values = c("#4393C3", "#D6604D")) +
      scale_y_continuous(parse(text = paste0("Mother~", out_label)), 
                         limits = c(lower_limit, upper_limit)) +
      scale_x_continuous(parse(text = paste0("Secondary~Children~", out_label)), 
                         limits = c(lower_limit, upper_limit)) +
      facet_wrap(~virus) +
      tmp_theme
    
    pl2 = exposure_data_summary %>%
      ungroup() %>%
      mutate(
        obs_infected = factor(obs_infected),
        idpar = fct_recode(fct_rev(idpar), "Secondary\nChildren" = "S", "Mother" = "M", "Household\nSum" = "HH")
        ) %>%
      ggplot(aes_string(x = "idpar", y = s, colour = "obs_infected")) +
      geom_boxplot() +
      scale_colour_manual(values = c("#4393C3", "#D6604D")) +
      geom_point(position = position_dodge(width = 0.75)) +
      scale_y_continuous(parse(text = out_label), limits = c(lower_limit, upper_limit)) +
      xlab(parse(text = c(""))) +
      facet_wrap(~virus) +
      tmp_theme +theme(axis.text.x = element_text(size = 7))
    
    plot_grid(pl1, pl2, nrow = 1, labels = LETTERS[c(pl_lab$letter_code, pl_lab$letter_code + 3)])
    
})

plot_grid(plot_grid(plotlist = pl_exposure, nrow = 3), trans_legend, nrow = 2, rel_heights = c(11, 1))

Past versions of exposure-profile-1.png

exposure_data_summary %>%
  subset(idpar != "HH") %>%
  group_by(virus, FamilyID, idpar) %>%
  gather(outcome, value, pct_pos, mean, maximum) %>%
  select(virus, FamilyID, idpar, outcome, value) %>%
  spread(idpar, value) %>%
  ungroup() %>%
  dplyr::mutate(
    fid = factor(FamilyID)
  ) %>%
  group_by(virus, outcome) %>%
  nest() %>%
  mutate(
    cor_test = map(data, ~cor.test(.x$M, .x$S), method = "Spearmen"),
    cor_res = map(cor_test, broom::tidy)
  ) %>%
  unnest(cor_res) %>%
  mutate(
    Endpoint = factor(outcome, levels = c("pct_pos", "mean", "maximum"),
                      labels = c("Pct. Positive", "Mean VL", "Maximum VL")),
    Estimate = stat_paste(estimate, conf.low, conf.high, digits = 2)
  ) %>%
  select(virus, Endpoint, Estimate, p.value) %>%
  arrange(virus, Endpoint) %>%
  kable(digits = 3) %>%
  kable_styling(full_width = F) %>%
  collapse_rows(1)

exposure_data_tests = 
  exposure_data_summary %>% 
  gather(stat, estimate, mean, maximum, pct_pos) %>%
  group_by(virus, stat, idpar) %>%
  mutate(obs_infected = factor(obs_infected)) %>%
  summarize(
    wilcox_pvalue =  coin::pvalue(coin::wilcox_test(estimate ~ obs_infected, distribution = "exact"))
  )


overall_summary = exposure_data_summary %>% 
  gather(stat, estimate, mean, maximum, pct_pos) %>%
  group_by(virus, stat, idpar, obs_infected) %>%
  summarize(
    median = median(estimate)
  ) %>%
  mutate(obs_infected = recode_factor(obs_infected, `1` = "Transmission", `0` = "No transmission", .ordered = T)) %>%
  spread(obs_infected, median) %>%
  mutate(Difference = `Transmission`-`No transmission`) %>%
  left_join(exposure_data_tests)

Joining, by = c("virus", "stat", "idpar")

overall_summary %>%
  kable(digits = 3) %>%
  kable_styling(full_width = F) %>%
  collapse_rows(columns = 1:2)

Household composition

Household sum composition was determined and reported by taking the mean proportion over measurements within a household. The summary across the households uses median and IQR to match the box plot statistics.

hh_summary = exposure_data %>%
  filter(HH > 0) %>%
  mutate(   
    S_pctHH = if_else(HH == 0, 0, 100 * (10^S/10^HH)),
    M_pctHH = if_else(HH == 0, 0, 100 * (10^M/10^HH))
    ) %>%
  group_by(virus, FamilyID, obs_infected) %>%
  summarise_at(vars(S_pctHH, M_pctHH), funs(mean = mean, median = median))

Warning: funs() is soft deprecated as of dplyr 0.8.0
please use list() instead

  # Before:
  funs(name = f(.))

  # After: 
  list(name = ~ f(.))
This warning is displayed once per session.

hh_summary %>%
  ungroup() %>%
  select(-obs_infected) %>%
  gather(stat, est, -virus, -FamilyID) %>%
  group_by(virus, stat) %>%
  summarize_if(is.double, funs(mean, median, IQR = IQR_range_str, range = range_str)) %>%
  mutate_if(is.double, round, digits = 2) %>%
  filter(str_detect(stat, "mean")) %>%
  mutate(
    stat = substr(stat, 1, 1)
  ) %>%
  rename(`Household member` = stat) %>%
  rename_at(vars(-IQR), funs(str_to_title)) %>%
  kable(caption = "Percent household composition") %>% kable_styling(full_width = F)

`mutate_if()` ignored the following grouping variables:
Column `virus`

hh_summary %>%
  ggplot(aes(x = virus, y = S_pctHH_mean)) +
  geom_boxplot() +
  geom_point(aes(colour = factor(obs_infected))) +
  xlab("") +
  scale_colour_manual("", values = infection_labels$colours, breaks = infection_labels$breaks,
                      labels = infection_labels$labels) +
  ylab("Percent of household shedding attributable to secondary children") 

Past versions of hh-composition-1.png

Interpolated exposure

In the sensitivity analysis, we’d like to assess the interpolation. All of the summary statistical analysis is limited to the observed exposures. Assuming the interpolation was unbiased, the precision in tests could still be artificially inflated without some correction. Because linear interpolation is used, the mean and maximum exposure estimates should be largely unaffected. However, the percent positive estimate may not be precise without a larger sample size.

exposure_interpolated_summary = exposure_data_long %>%
  mutate(
    pos_count = count > 0
    ) %>%
  group_by(virus, FamilyID, obs_infected, idpar) %>%
  mutate(
    interpolated_pct = 100*mean(interpolated),
    total_pos = sum(pos_count),
    pct_pos = 100 * mean(pos_count)
    ) %>%
  group_by(pct_pos, total_pos, interpolated_pct, add = T) %>%
  summarise_at(vars(count), funs(N = n(), mean = mean, median = median, maximum = max))

exposure_interpolated_summary %>%
  ungroup() %>%
  mutate(
    idpar = fct_recode(fct_rev(idpar), 
                       "Secondary\nChildren" = "S", "Mother" = "M", "Household\nSum" = "HH")
    ) %>%
  ggplot(aes(x = idpar, y = interpolated_pct)) +
  geom_boxplot() +
  geom_point() +
  xlab("") +
  ylab("Percent of weekly exposures interpolated") +
  facet_wrap(~virus)

Past versions of interpolate-plot-1.png

exposure_interpolated_summary %>% 
  ungroup() %>%
  subset(idpar != "HH") %>%
  gather(stat, estimate, mean, maximum, pct_pos) %>%
  mutate(
    idpar = fct_recode(fct_rev(idpar), 
                       "Secondary\nChildren" = "S", "Mother" = "M"),
    stat = fct_recode(fct_rev(stat), 
                       "% positive (weekly)" = "pct_pos", "Mean VL" = "mean", "Maximum VL" = "maximum"),
    ) %>%
  ggplot(aes(x = interpolated_pct, y = estimate, colour = factor(obs_infected))) +
  geom_point() +
  ylab("") +
  xlab("Percent of weekly exposures interpolated") +
  facet_grid(stat ~ virus+idpar, scales = "free_y", switch = "y") +
  scale_colour_manual("", values = infection_labels$colours, breaks = infection_labels$breaks,
                      labels = infection_labels$labels) +
  theme(strip.placement = "outside")

Past versions of interpolate-exposure-1.png

Session information

sessionInfo()

R version 3.6.0 (2019-04-26)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS Mojave 10.14.5

Matrix products: default
BLAS:   /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] cowplot_0.9.4     kableExtra_1.1.0  conflicted_1.0.3 
 [4] forcats_0.4.0     stringr_1.4.0     dplyr_0.8.1      
 [7] purrr_0.3.2       readr_1.3.1       tidyr_0.8.3      
[10] tibble_2.1.3      ggplot2_3.1.1     tidyverse_1.2.1  
[13] survival_2.44-1.1

loaded via a namespace (and not attached):
 [1] httr_1.4.0         jsonlite_1.6       viridisLite_0.3.0 
 [4] splines_3.6.0      modelr_0.1.4       assertthat_0.2.1  
 [7] highr_0.8          stats4_3.6.0       selectr_0.4-1     
[10] coin_1.3-0         cellranger_1.1.0   yaml_2.2.0        
[13] pillar_1.4.1       backports_1.1.4    lattice_0.20-38   
[16] glue_1.3.1         digest_0.6.19      rvest_0.3.4       
[19] colorspace_1.4-1   sandwich_2.5-1     htmltools_0.3.6   
[22] Matrix_1.2-17      plyr_1.8.4         pkgconfig_2.0.2   
[25] broom_0.5.2        haven_2.1.0        mvtnorm_1.0-11    
[28] scales_1.0.0       webshot_0.5.1      whisker_0.3-2     
[31] git2r_0.25.2       generics_0.0.2     TH.data_1.0-10    
[34] withr_2.1.2        lazyeval_0.2.2     cli_1.1.0         
[37] magrittr_1.5       crayon_1.3.4       readxl_1.3.1      
[40] memoise_1.1.0      evaluate_0.14      fs_1.3.1          
[43] nlme_3.1-140       MASS_7.3-51.4      xml2_1.2.0        
[46] tools_3.6.0        hms_0.4.2          matrixStats_0.54.0
[49] multcomp_1.4-10    munsell_0.5.0      compiler_3.6.0    
[52] rlang_0.4.0        grid_3.6.0         rstudioapi_0.10   
[55] labeling_0.3       rmarkdown_1.13     codetools_0.2-16  
[58] gtable_0.3.0       reshape2_1.4.3     R6_2.4.0          
[61] zoo_1.8-6          lubridate_1.7.4    knitr_1.23        
[64] workflowr_1.3.0    libcoin_1.0-4      rprojroot_1.3-2   
[67] modeltools_0.2-22  stringi_1.4.3      parallel_3.6.0    
[70] Rcpp_1.0.1         tidyselect_0.2.5   xfun_0.7