class: center, middle, inverse, title-slide # serosurvey:<br/>Serological Surveys and<br/>Prevalence Estimation<br/>Under Misclassification ### Andree Valle-Campos (<span class="citation">@avallecam</span>) ### Centro Nacional de Epidemiología, CDC Perú --- class: center, middle <style type="text/css"> pre { max-width: 200%; overflow-x: scroll; } </style> # slides:<br/>[bit.ly/serosurvey-user21](https://bit.ly/serosurvey-user21) --- # Introduction -- - Population-based serological surveys __generate__ seroprevalence estimates and __quantify__ how many people have been infected by a certain pathogen. -- - These studies need to estimate prevalence considering __misclassification__ due to imperfect diagnostic tests with sensitivity and specificity __known with certainty__. -- - However, during __first months__ of a novel pathogen outbreak, __performance__ of diagnostic tests are __“unknown” or “uncertain”__ given limited validation studies. -- - This situation opened the opportunity for the development of newer __ad-hoc methods__. Most of them used R, but not stored in a __R package__. -- - The __increased demand__ for prevalence estimates and the absence of standardized procedures favored the usage of __non-reproducible workflows__. -- - We created the __[serosurvey](https://avallecam.github.io/serosurvey/) R package__ to gather known and recently developed __functions__ and build __workflow__ templates. --- class: chapter-slide # The `serosurvey` functions --- class: scrollable-slide <!-- `survey` --> # Estimate single prevalence -- From a [`srvyr`](http://gdfe.co/srvyr/) [`survey`](https://r-survey.r-forge.r-project.org/survey/) design object, __`serosvy_proportion`__ estimates: + weighted prevalence (`prop`), + total population (`total`), + raw proportion (`raw_prop`), + coefficient of variability (`cv`), + design effect (`deff`) -- ```r *serosvy_proportion(design = design, denominator = covariate_01, numerator = outcome_one) ``` ``` ## # A tibble: 6 x 23 ## denominator denominator_lev~ numerator numerator_level prop prop_low prop_upp ## <chr> <fct> <chr> <fct> <dbl> <dbl> <dbl> ## 1 covariate_~ E outcome_~ No 0.211 0.130 0.325 ## 2 covariate_~ E outcome_~ Yes 0.789 0.675 0.870 ## 3 covariate_~ H outcome_~ No 0.852 0.564 0.962 ## 4 covariate_~ H outcome_~ Yes 0.148 0.0377 0.436 ## 5 covariate_~ M outcome_~ No 0.552 0.224 0.840 ## 6 covariate_~ M outcome_~ Yes 0.448 0.160 0.776 ## # ... with 16 more variables: prop_cv <dbl>, prop_se <dbl>, total <dbl>, ## # total_low <dbl>, total_upp <dbl>, total_cv <dbl>, total_se <dbl>, ## # total_deff <dbl>, total_den <dbl>, total_den_low <dbl>, ## # total_den_upp <dbl>, raw_num <int>, raw_den <int>, raw_prop <dbl>, ## # raw_prop_low <dbl>, raw_prop_upp <dbl> ``` <!-- `serology` --> <!-- misclassification --> <!-- - We gather __one frequentist approach__ (Rogan and Gladen, 1978), --> <!-- available in Github repos (Azman et al. 2020; Takahashi et al. 2020). --> <!-- Check [here](https://avallecam.github.io/serosurvey/reference/index.html). --> --- # Known test performance -- - We provide functions with __tidy outputs__ for frequentist and bayesian approaches available in GitHub repos to solve the inherited __misclassification__. -- - One function adapted the __bayesian approach__ developed by Larremore, et al. 2021 for [known](https://github.com/LarremoreLab/covid_serological_sampling/blob/master/codebase/seroprevalence.R) test performance -- ```r *serosvy_known_sample_posterior( #in population positive_number_test = 321, total_number_test = 321+1234, # known performance sensitivity = 0.93, specificity = 0.975) ``` ``` ## # A tibble: 1 x 2 ## posterior summary ## <list> <list> ## 1 <tibble [10,000 x 3]> <tibble [1 x 5]> ``` ``` ## # A tibble: 1 x 5 ## skim_variable numeric.p05 numeric.mean numeric.p50 numeric.p95 ## <chr> <dbl> <dbl> <dbl> <dbl> ## 1 r1 0.183 0.201 0.201 0.219 ``` <img src="serosurvey-user21_files/figure-html/unnamed-chunk-9-1.png" title="One ggplot panel with density as vertical y-axis. Histogram and density plot of the Posterior Distribution of the Prevalence accounting for misclassification from a diagnostic test with known performance. Additionally, it includes vertical lines at the mean, 5% and 95% percentiles." alt="One ggplot panel with density as vertical y-axis. Histogram and density plot of the Posterior Distribution of the Prevalence accounting for misclassification from a diagnostic test with known performance. Additionally, it includes vertical lines at the mean, 5% and 95% percentiles." /> --- # Unknown test performance -- - And a second one did the same with other __bayesian approach__ developed by Larremore et al. 2020 for [unknown](https://github.com/LarremoreLab/bayesian-joint-prev-se-sp/blob/master/singleSERO_uncertainTEST.R) test performance. -- ```r *serosvy_unknown_sample_posterior( #in population positive_number_test = 321, total_number_test = 321+1234, # in lab (local validation study) true_positive = 670, true_negative = 640, false_positive = 202, false_negative = 74) ``` ``` ## # A tibble: 1 x 4 ## posterior summary performance tidy ## <list> <list> <list> <list> ## 1 <tibble [10,000 x 3]> <tibble [1 x 5]> <tibble [2 x 5]> <tibble [3 x 5]> ``` ``` ## # A tibble: 3 x 5 ## skim_variable numeric.p05 numeric.mean numeric.p50 numeric.p95 ## <chr> <dbl> <dbl> <dbl> <dbl> ## 1 r 0.000674 0.0100 0.00757 0.0274 ## 2 se 0.881 0.899 0.900 0.917 ## 3 sp 0.771 0.786 0.786 0.801 ``` <img src="serosurvey-user21_files/figure-html/unnamed-chunk-13-1.png" title="Three ggplot facet panels with the same vertical y-axis of density. Histogram and density plots of the Posterior Distribution of the (i) r: Prevalence accounting for misclassification from a diagnostic test with unknown performance; (ii) se: Sensitivity; (iii) sp: Specificity." alt="Three ggplot facet panels with the same vertical y-axis of density. Histogram and density plots of the Posterior Distribution of the (i) r: Prevalence accounting for misclassification from a diagnostic test with unknown performance; (ii) se: Sensitivity; (iii) sp: Specificity." /> --- class: chapter-slide # The `serosurvey` workflow --- # Multiple prevalences -- - We built a [workflow](https://avallecam.github.io/serosurvey/articles/howto-reprex.html) to __estimate multiple prevalences__ + suitable for different sets of numerators and denominators, + in one single pipe operation -- ```r # set 01 of denominator-numerator *expand_grid( design=list(design), denominator=c("covariate_01","covariate_02"), # covariates numerator=c("outcome_one","outcome_two") # outcomes ) %>% # set 02 of denominator-numerator (e.g. within main outcome) * union_all( expand_grid( design=list(design), denominator=c("outcome_one","outcome_two"), # outcomes numerator=c("covariate_02") # covariates )) %>% # create symbols (to be read as arguments) mutate( * denominator=map(denominator,dplyr::sym), numerator=map(numerator,dplyr::sym) ) %>% # estimate prevalence mutate(output=pmap(.l = select(.,design,denominator,numerator), * .f = serosvy_proportion)) %>% # show the outcome select(-design,-denominator,-numerator) %>% unnest(cols = c(output)) %>% print(n=Inf) ``` ``` ## # A tibble: 25 x 23 ## denominator denominator_level numerator numerator_level prop prop_low ## <chr> <fct> <chr> <fct> <dbl> <dbl> ## 1 covariate_01 E outcome_one No 0.211 0.130 ## 2 covariate_01 E outcome_one Yes 0.789 0.675 ## 3 covariate_01 H outcome_one No 0.852 0.564 ## 4 covariate_01 H outcome_one Yes 0.148 0.0377 ## 5 covariate_01 M outcome_one No 0.552 0.224 ## 6 covariate_01 M outcome_one Yes 0.448 0.160 ## 7 covariate_01 E outcome_two (-0.1,50] 0.182 0.0499 ## 8 covariate_01 E outcome_two (50,100] 0.818 0.515 ## 9 covariate_01 H outcome_two (-0.1,50] 0.0769 0.00876 ## 10 covariate_01 H outcome_two (50,100] 0.923 0.560 ## 11 covariate_01 M outcome_two (50,100] 1.00 1.00 ## 12 covariate_02 No outcome_one No 1.00 1.00 ## 13 covariate_02 Yes outcome_one No 0.0334 0.00884 ## 14 covariate_02 Yes outcome_one Yes 0.967 0.882 ## 15 covariate_02 No outcome_two (-0.1,50] 0.218 0.0670 ## 16 covariate_02 No outcome_two (50,100] 0.782 0.479 ## 17 covariate_02 Yes outcome_two (-0.1,50] 0.0914 0.0214 ## 18 covariate_02 Yes outcome_two (50,100] 0.909 0.684 ## 19 outcome_one No covariate_02 No 0.939 0.778 ## 20 outcome_one No covariate_02 Yes 0.0615 0.0148 ## 21 outcome_one Yes covariate_02 Yes 1.00 1.00 ## 22 outcome_two (-0.1,50] covariate_02 No 0.549 0.294 ## 23 outcome_two (-0.1,50] covariate_02 Yes 0.451 0.219 ## 24 outcome_two (50,100] covariate_02 No 0.305 0.188 ## 25 outcome_two (50,100] covariate_02 Yes 0.695 0.546 ## # ... with 17 more variables: prop_upp <dbl>, prop_cv <dbl>, prop_se <dbl>, ## # total <dbl>, total_low <dbl>, total_upp <dbl>, total_cv <dbl>, ## # total_se <dbl>, total_deff <dbl>, total_den <dbl>, total_den_low <dbl>, ## # total_den_upp <dbl>, raw_num <int>, raw_den <int>, raw_prop <dbl>, ## # raw_prop_low <dbl>, raw_prop_upp <dbl> ``` --- # Add test performance -- - Using [`purrr`](https://purrr.tidyverse.org/) and [`furrr`](https://davisvaughan.github.io/furrr/) within a `for-loop` to efficiently iterate and parallelize the __update of prevalence posterior distributions__ -- - We apply the procedure to __each__ point estimate and confidence interval -- ```r # 56-60sec per covariate # 4GB RAM # parallelized in 8 cores using purrr y furrr *plan(multisession, workers = availableCores()) tic() out <- tibble() *for (i in 1:nrow(outcome_01_adj_pre)) { out <- outcome_01_adj_pre %>% slice(i) %>% # dot estimate mutate(adj_dot_unk=future_pmap( .l = select(., positive_number_test=total_round, total_number_test=total_den_round), * .f = possibly(serosvy_unknown_sample_posterior,otherwise = NA_real_), true_positive = true_positive, false_negative = false_negative, false_positive = false_positive, true_negative = true_negative)) %>% * serosvy_extract_posterior(variable = adj_dot_unk) %>% * # low confidence interval mutate(adj_low_unk=future_pmap( .l = select(., positive_number_test=total_low_round, total_number_test=total_den_low_round), .f = possibly(serosvy_unknown_sample_posterior,otherwise = NA_real_), true_positive = true_positive, false_negative = false_negative, false_positive = false_positive, true_negative = true_negative)) %>% serosvy_extract_posterior(variable = adj_low_unk) %>% * # upper confidence interval mutate(adj_upp_unk=future_pmap( .l = select(., positive_number_test=total_upp_round, total_number_test=total_den_upp_round), .f = possibly(serosvy_unknown_sample_posterior,otherwise = NA_real_), true_positive = true_positive, false_negative = false_negative, false_positive = false_positive, true_negative = true_negative)) %>% serosvy_extract_posterior(variable = adj_upp_unk) %>% # union all outputs union_all(out) # print the progress! (recommended) * out %>% print() } toc() *outcome_01_adj <- out %>% mutate(rowname=as.numeric(rowname)) %>% arrange(rowname) %>% select(-rowname) ``` --- # Take home messages -- - The `serosurvey` R package facilitated the estimation of prevalence in the __current epidemic__. -- - The `serosurvey` R package wanted to reduce the __reproducibility gap__ of serological survey analysis using diagnostic tests with unknown performance within a free software environment. -- - The usage of an __R package plus a workflow__ may improve our preparedness and response for next epidemics. -- - Specially in our capacity to __generate comparable results__ from studies with similar epidemiological designs and statistical analysis plans. --- class: center, middle # Thanks for your attention! <br/> Happy to chat with you! ## <br/>Andree Valle Campos<br/><br/> <svg viewBox="0 0 512 512" style="height:1em;position:relative;display:inline-block;top:.1em;" xmlns="http://www.w3.org/2000/svg"> <path d="M459.37 151.716c.325 4.548.325 9.097.325 13.645 0 138.72-105.583 298.558-298.558 298.558-59.452 0-114.68-17.219-161.137-47.106 8.447.974 16.568 1.299 25.34 1.299 49.055 0 94.213-16.568 130.274-44.832-46.132-.975-84.792-31.188-98.112-72.772 6.498.974 12.995 1.624 19.818 1.624 9.421 0 18.843-1.3 27.614-3.573-48.081-9.747-84.143-51.98-84.143-102.985v-1.299c13.969 7.797 30.214 12.67 47.431 13.319-28.264-18.843-46.781-51.005-46.781-87.391 0-19.492 5.197-37.36 14.294-52.954 51.655 63.675 129.3 105.258 216.365 109.807-1.624-7.797-2.599-15.918-2.599-24.04 0-57.828 46.782-104.934 104.934-104.934 30.213 0 57.502 12.67 76.67 33.137 23.715-4.548 46.456-13.32 66.599-25.34-7.798 24.366-24.366 44.833-46.132 57.827 21.117-2.273 41.584-8.122 60.426-16.243-14.292 20.791-32.161 39.308-52.628 54.253z"></path></svg> <svg viewBox="0 0 496 512" style="height:1em;position:relative;display:inline-block;top:.1em;" xmlns="http://www.w3.org/2000/svg"> <path d="M165.9 397.4c0 2-2.3 3.6-5.2 3.6-3.3.3-5.6-1.3-5.6-3.6 0-2 2.3-3.6 5.2-3.6 3-.3 5.6 1.3 5.6 3.6zm-31.1-4.5c-.7 2 1.3 4.3 4.3 4.9 2.6 1 5.6 0 6.2-2s-1.3-4.3-4.3-5.2c-2.6-.7-5.5.3-6.2 2.3zm44.2-1.7c-2.9.7-4.9 2.6-4.6 4.9.3 2 2.9 3.3 5.9 2.6 2.9-.7 4.9-2.6 4.6-4.6-.3-1.9-3-3.2-5.9-2.9zM244.8 8C106.1 8 0 113.3 0 252c0 110.9 69.8 205.8 169.5 239.2 12.8 2.3 17.3-5.6 17.3-12.1 0-6.2-.3-40.4-.3-61.4 0 0-70 15-84.7-29.8 0 0-11.4-29.1-27.8-36.6 0 0-22.9-15.7 1.6-15.4 0 0 24.9 2 38.6 25.8 21.9 38.6 58.6 27.5 72.9 20.9 2.3-16 8.8-27.1 16-33.7-55.9-6.2-112.3-14.3-112.3-110.5 0-27.5 7.6-41.3 23.6-58.9-2.6-6.5-11.1-33.3 2.6-67.9 20.9-6.5 69 27 69 27 20-5.6 41.5-8.5 62.8-8.5s42.8 2.9 62.8 8.5c0 0 48.1-33.6 69-27 13.7 34.7 5.2 61.4 2.6 67.9 16 17.7 25.8 31.5 25.8 58.9 0 96.5-58.9 104.2-114.8 110.5 9.2 7.9 17 22.9 17 46.4 0 33.7-.3 75.4-.3 83.6 0 6.5 4.6 14.4 17.3 12.1C428.2 457.8 496 362.9 496 252 496 113.3 383.5 8 244.8 8zM97.2 352.9c-1.3 1-1 3.3.7 5.2 1.6 1.6 3.9 2.3 5.2 1 1.3-1 1-3.3-.7-5.2-1.6-1.6-3.9-2.3-5.2-1zm-10.8-8.1c-.7 1.3.3 2.9 2.3 3.9 1.6 1 3.6.7 4.3-.7.7-1.3-.3-2.9-2.3-3.9-2-.6-3.6-.3-4.3.7zm32.4 35.6c-1.6 1.3-1 4.3 1.3 6.2 2.3 2.3 5.2 2.6 6.5 1 1.3-1.3.7-4.3-1.3-6.2-2.2-2.3-5.2-2.6-6.5-1zm-11.4-14.7c-1.6 1-1.6 3.6 0 5.9 1.6 2.3 4.3 3.3 5.6 2.3 1.6-1.3 1.6-3.9 0-6.2-1.4-2.3-4-3.3-5.6-2z"></path></svg> __`@avallecam`__ <br/> <svg viewBox="0 0 576 512" style="height:1em;position:relative;display:inline-block;top:.1em;" xmlns="http://www.w3.org/2000/svg"> <path d="M567.938 243.908L462.25 85.374A48.003 48.003 0 0 0 422.311 64H153.689a48 48 0 0 0-39.938 21.374L8.062 243.908A47.994 47.994 0 0 0 0 270.533V400c0 26.51 21.49 48 48 48h480c26.51 0 48-21.49 48-48V270.533a47.994 47.994 0 0 0-8.062-26.625zM162.252 128h251.497l85.333 128H376l-32 64H232l-32-64H76.918l85.334-128z"></path></svg> __avallecam@gmail.com__ --- class: chapter-slide # Extra slides --- # References (1/3) [1] Takahashi, S., et al. 2020. Are SARS-CoV-2 seroprevalence estimates biased?. _The Journal of Infectious Diseases_. https://doi.org/10.1093/infdis/jiaa523 [2] Flor, M., et al. 2020. “Comparison of Bayesian and Frequentist Methods for Prevalence Estimation Under Misclassification.” BMC Public Health 20 (1). https://doi.org/10.1186/s12889-020-09177-4. [3] Kritsotakis, EI., 2020. “On the Importance of Population-Based Serological Surveys of SARS-CoV-2 Without Overlooking Their Inherent Uncertainties.” _Public Health in Practice_ 1 (November): 100013. [10.1016/j.puhip.2020.100013](https://doi.org/10.1016/j.puhip.2020.100013). --- # References (2/3) [4] Larremore, DB., Bailey K. Fosdick, Sam Zhang, and Yonatan H. Grad. 2020. “Jointly Modeling Prevalence, Sensitivity and Specificity for Optimal Sample Allocation.” _bioRxiv_, May. https://doi.org/10.1101/2020.05.23.112649. [5] Gelman, Andrew, and Bob Carpenter. 2020. “Bayesian Analysis of Tests with Unknown Specificity and Sensitivity.” _Journal of the Royal Statistical Society: Series C (Applied Statistics)_, August. https://doi.org/10.1111/rssc.12435. [6] Andree Valle Campos. 2020. avallecam/serosurvey: Serological Survey Analysis For Prevalence Estimation Under Misclassification (Version v1.0.1). Zenodo. http://doi.org/10.5281/zenodo.4118710. https://avallecam.github.io/serosurvey --- # References (3/3) [7] Rogan, Walter J., and Beth Gladen. 1978. “Estimating Prevalence from the Results of A Screening Test.” _American Journal of Epidemiology_ 107 (1): 71–76. https://doi.org/10.1093/oxfordjournals.aje.a112510. [8] Larremore, DB, et al. 2021. “Estimating SARS-CoV-2 Seroprevalence and Epidemiological Parameters with Uncertainty from Serological Surveys.” _eLife_, March. https://doi.org/10.7554/eLife.64206. [9] Larremore, DB., et al. 2020. “Jointly Modeling Prevalence, Sensitivity and Specificity for Optimal Sample Allocation.” _bioRxiv_, May. https://doi.org/10.1101/2020.05.23.112649. --- # Key Results -- - For the __survey__ side, we provide functions to calculate single prevalences and generate __tidy outputs__ for a hierarchical bayesian approach that incorporates the __unknown__ test performance -- - For the __serology__ side, we gather functions of methods available in GitHub repos to solve the inherited misclassification in population-based serological surveys. -- - We created a __workflow__ that simulates the usage of an __imperfect test__ to measure an outcome within a __multi-stage sampling__ survey design and incorporate the test uncertainty into the sampling design uncertainty. -- - The `serosurvey` R package was used in the analysis code of the Lima seroprevalence study. -- [10] Reyes-Vega, Mary F., et al. "SARS-CoV-2 prevalence associated to low socioeconomic status and overcrowding in an LMIC megacity: A population-based seroepidemiological survey in Lima, Peru." _EClinicalMedicine_ 34 (2021): 100801. https://doi.org/10.1016/j.eclinm.2021.100801 --- # Next steps -- - Adjust to the __rOpenSci__ packaging guidelines -- - Submit package to __CRAN__ -- - Integrate it to the R Epidemic Consortium __RECON__ -- - Welcome new functions and workflows to the package via __GitHub__! -- .footnote[ [1] https://devguide.ropensci.org/building.html [2] https://www.repidemicsconsortium.org/resources/guidelines/ [3] https://twitter.com/avallecam/status/1312670190694535169 ]