diff --git a/.Rbuildignore b/.Rbuildignore index f962eb22..958b4f56 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -23,3 +23,10 @@ allpopsamples_hlye.csv$ ^CRAN-SUBMISSION$ ^README\.qmd$ ^codecov\.yml$ +^_quarto\.yml$ +^\.lintr$ +^vignettes$ +^vignettes/\.quarto$ +^vignettes/methodology\.qmd$ +^\.quarto$ +^man/df_to_array\.Rd$ diff --git a/.gitignore b/.gitignore index 02f621ff..f2f0195f 100644 --- a/.gitignore +++ b/.gitignore @@ -15,3 +15,5 @@ serocalculator*.tar.gz serocalculator*.tgz README.html README_files + +/.quarto/ diff --git a/.lintr b/.lintr new file mode 100644 index 00000000..3d4eb22e --- /dev/null +++ b/.lintr @@ -0,0 +1,2 @@ +linters: linters_with_defaults() # see vignette("lintr") +encoding: "UTF-8" diff --git a/DESCRIPTION b/DESCRIPTION index e9461298..669cd35a 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Type: Package Package: serocalculator Title: Estimating Infection Rates from Serological Data -Version: 1.2.0.9009 +Version: 1.2.0.9014 Authors@R: c( person("Peter", "Teunis", , "p.teunis@emory.edu", role = c("aut", "cph"), comment = "Author of the method and original code."), @@ -41,22 +41,31 @@ Imports: utils Suggests: bookdown, + devtag (>= 0.0.0.9000), DT, fs, ggbeeswarm, knitr, + pak, parallel, readr, + quarto, rmarkdown, spelling, + ssdtools (>= 1.0.6.9016), testthat (>= 3.0.0), - tidyverse + tidyverse, + qrcode LinkingTo: Rcpp Config/testthat/edition: 3 +Config/Needs/build: moodymudskipper/devtag Encoding: UTF-8 Language: en-US LazyData: true NeedsCompilation: no -Roxygen: list(markdown = TRUE) +Remotes: + bcgov/ssdtools, + moodymudskipper/devtag +Roxygen: list(markdown = TRUE, roclets = c("collate", "rd", "namespace", "devtag::dev_roclet")) RoxygenNote: 7.3.2 diff --git a/NAMESPACE b/NAMESPACE index b94ff416..433fd515 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -6,23 +6,10 @@ S3method(autoplot,pop_data) S3method(autoplot,seroincidence) S3method(autoplot,seroincidence.by) S3method(autoplot,summary.seroincidence.by) -S3method(get_age,pop_data) -S3method(get_age_var,pop_data) -S3method(get_biomarker_levels,pop_data) -S3method(get_biomarker_names,pop_data) -S3method(get_biomarker_names_var,pop_data) -S3method(get_id,pop_data) -S3method(get_id_var,pop_data) -S3method(get_value,pop_data) -S3method(get_value_var,pop_data) S3method(print,seroincidence) S3method(print,seroincidence.by) S3method(print,summary.pop_data) S3method(print,summary.seroincidence.by) -S3method(set_age,pop_data) -S3method(set_biomarker_var,pop_data) -S3method(set_id,pop_data) -S3method(set_value,pop_data) S3method(strata,seroincidence.by) S3method(summary,pop_data) S3method(summary,seroincidence) @@ -30,6 +17,7 @@ S3method(summary,seroincidence.by) export(as_curve_params) export(as_noise_params) export(as_pop_data) +export(autoplot) export(check_pop_data) export(est.incidence) export(est.incidence.by) @@ -39,7 +27,7 @@ export(fdev) export(getAdditionalData) export(get_additional_data) export(graph.curve.params) -export(graph.loglik) +export(graph_loglik) export(llik) export(load_curve_params) export(load_noise_params) diff --git a/NEWS.md b/NEWS.md index 4dbeb135..25032fbf 100644 --- a/NEWS.md +++ b/NEWS.md @@ -7,6 +7,14 @@ ## New features +* Updated `simulate_xsectionalData.Rmd()` (linting, removing deprecated functions). + +* Added default value for `antigen_isos` argument in `log_likelihood()` (#286) + +* Updated enteric fever example article with upgraded code and visualizations (#290) + +* Added `Methodology` vignette (#284, #302) + * Added template for reporting Issues (from `usethis::use_tidy_issue_template()`) (#270) @@ -15,11 +23,23 @@ ## Developer-facing changes -* Updated GitHub Action files and reformatted `DESCRIPTION` (#268) +* initialized [`lintr`](https://lintr.r-lib.org/) with `lintr::use_lint()` (#278) +* created unit test for `df_to_array()` (#276) + +* fixed `dplyr::select()` deprecation warning in `df_to_array()` (#276) + +* Added `devtag` to package (using `devtag::use_devtag()`) (#292) +* Added `@dev` tag to `?df_to_array()` (#292) + +* Generalized `get_()` and `set_()` methods to be general-purpose +(no S3 class-specific methods needed yet) (#274). + +* Updated GitHub Action files and reformatted `DESCRIPTION` (#268) * Added `.gitattributes` file () copied from +* Added QR code to `README.qmd` * Added additional automated checks through [GitHub actions](https://docs.github.com/en/actions), including: diff --git a/R/as_curve_params.R b/R/as_curve_params.R index 869a7370..6dac7101 100644 --- a/R/as_curve_params.R +++ b/R/as_curve_params.R @@ -1,8 +1,10 @@ #' Load antibody decay curve parameter #' #' @param data a [data.frame()] or [tibble::tbl_df] -#' @param antigen_isos [character()] vector of antigen isotypes to be used in analyses -#' @returns a `curve_data` object (a [tibble::tbl_df] with extra attribute `antigen_isos`) +#' @param antigen_isos a [character()] vector of antigen isotypes +#' to be used in analyses +#' @returns a `curve_data` object +#' (a [tibble::tbl_df] with extra attribute `antigen_isos`) #' @export #' @examples #' library(magrittr) @@ -14,8 +16,7 @@ #' print(curve_data) as_curve_params <- function(data, antigen_isos = NULL) { - if(!is.data.frame(data)) - { + if (!is.data.frame(data)) { cli::cli_abort( class = "not data.frame", message = c( @@ -31,16 +32,18 @@ as_curve_params <- function(data, antigen_isos = NULL) { data %>% tibble::as_tibble() + # check if object has expected columns: + # define curve columns curve_cols <- c("antigen_iso", "y0", "y1", "t1", "alpha", "r") - # check if object is curve (with columns) - if (!all(is.element(curve_cols, curve_data %>% names()))) { - # get columns from provided data - data_cols <- data %>% names() + # get columns from provided data + data_cols <- data %>% names() + + # get any missing column(s) + missing_cols <- setdiff(x = curve_cols, y = data_cols) - # get any missing column(s) - missing_cols <- setdiff(x = curve_cols, y = data_cols) + if (length(missing_cols) > 0) { cli::cli_abort( class = "not curve_params", @@ -67,5 +70,8 @@ as_curve_params <- function(data, antigen_isos = NULL) { # assign antigen attribute attr(curve_data, "antigen_isos") <- antigen_isos + curve_data <- curve_data %>% + set_biomarker_var(biomarker = "antigen_iso", standardize = FALSE) + return(curve_data) } diff --git a/R/class_attributes.R b/R/class_attributes.R new file mode 100644 index 00000000..19812672 --- /dev/null +++ b/R/class_attributes.R @@ -0,0 +1,191 @@ +get_age_var <- function(object, ...) { + age_var <- attributes(object)$age_var + return(age_var) +} + +get_age <- function(object, ...) { + age_var <- object %>% get_age_var() + age_data <- object %>% pull(age_var) + return(age_data) +} + +get_value_var <- function(object, ...) { + value_var <- attributes(object)$value_var + return(value_var) +} + +get_value <- function(object, ...) { + value_var_name <- object %>% get_value_var() + value_data <- object %>% pull(value_var_name) + return(value_data) +} + +get_id_var <- function(object, ...) { + id_var <- attributes(object)$id_var + return(id_var) +} + +get_id <- function(object, ...) { + id_var_name <- object %>% get_id_var() + id_data <- object %>% pull(id_var_name) + return(id_data) +} + +get_biomarker_levels <- function(object, ...) { + attr(object, "antigen_isos") +} + +get_biomarker_names_var <- function(object, ...) { + # get value attribute + biomarker_var <- attributes(object)[["biomarker_var"]] + + return(biomarker_var) +} + +get_biomarker_names <- function(object, ...) { + # get biomarker name data + biomarker_names_var <- get_biomarker_names_var(object) + biomarker_data <- object %>% pull(biomarker_names_var) + + return(biomarker_data) +} + + +set_age <- function(object, + age = "Age", + standardize = TRUE, + ...) { + # check if age column exists + if (age %in% colnames(object)) { + attr(object, "age_var") <- age + } else { + cli::cli_warn(class = "missing variable", + 'The specified `age` column "{age}" does not exist.') + + # search age variable from object + age_var <- + grep( + x = colnames(object), + value = TRUE, + pattern = age, + ignore.case = TRUE + ) + + if (length(age_var) == 1) { + attr(object, "age_var") <- age_var + + # create warning when using searched age instead of provided age + cli::cli_inform('Proceeding to use "{.var {age_var}}"') + } else if (length(age_var) == 0) { + cli::cli_abort("No similar column name was detected.") + } else if (length(age_var) > 1) { + cli::cli_warn("Multiple potential matches found: {.var {age_var}}") + cli::cli_warn("Using first match: {.var {age_var[1]}}") + attr(object, "age_var") <- age_var[1] + } else { + cli::cli_abort("{.code length(age_var)} = {.val {length(age_var)}}") + } + } + + if (standardize) { + object <- object %>% + rename(c("age" = attr(object, "age_var"))) + + # set age attribute + attr(object, "age_var") <- "age" + } + + return(object) +} + + +set_value <- function(object, + value = "result", + standardize = TRUE, + ...) { + # check if value column exists + if (value %in% colnames(object)) { + attr(object, "value_var") <- value + } else { + cli::cli_warn('The specified `value` column "{.var {value}}" + does not exist.') + + # search value variable from pop_data + value_var <- + grep( + x = colnames(object), + value = TRUE, + pattern = value, + ignore.case = TRUE + ) + + if (length(value_var) == 1) { + attr(object, "value_var") <- value_var + + # create warning when using searched age instead of provided age + cli::cli_inform('Proceeding to use "{.var {value_var}}"') + } else if (length(value_var) == 0) { + cli::cli_abort("No similar column name was detected.") + } else { + # i.e. if (length(value_var) > 1) + cli::cli_warn("Multiple potential matches found: {.var {value_var}}") + cli::cli_inform("Using first match: {.var {value_var[1]}}") + attr(object, "value_var") <- value_var[1] + } + } + + if (standardize) { + object <- object %>% + rename(c("value" = attr(object, "value_var"))) + + # set id attribute + attr(object, "value_var") <- "value" + } + + return(object) +} + +set_id <- function(object, + id = "index_id", + standardize = TRUE, + ...) { + # check if id column exists + if (id %in% colnames(object)) { + attr(object, "id_var") <- id + } else { + cli::cli_warn("The specified {.var id} column {.val {id}} does not exist.") + + # search id variable from object + id_var <- + grep( + x = colnames(object), + value = TRUE, + pattern = id, + ignore.case = TRUE + ) + + if (length(id_var) == 1) { + attr(object, "id_var") <- id_var + + # create warning when using searched id instead of provided id + cli::cli_inform('Proceeding to use "{id_var}"') + } else if (length(id_var) == 0) { + cli::cli_abort("No similar column name was detected.") + } else { + # if (length(id_var) > 1) + cli::cli_warn("Multiple potential matches found: {.var {id_var}}") + cli::cli_inform("Using first match: {.var {id_var[1]}}") + attr(object, "id_var") <- id_var[1] + } + } + + if (standardize) { + object <- object %>% + rename(c("id" = attr(object, "id_var"))) + + # set id attribute + attr(object, "id_var") <- "id" + } + + return(object) +} diff --git a/R/df_to_array.R b/R/df_to_array.R index 6c306240..98e5ff4b 100644 --- a/R/df_to_array.R +++ b/R/df_to_array.R @@ -8,7 +8,7 @@ #' #' @keywords internal #' -df.to.array <- function( +df.to.array <- function( # nolint: object_name_linter df, dim_var_names, value_var_name = "value") { @@ -18,10 +18,14 @@ df.to.array <- function( #' Convert a data.frame (or tibble) into a multidimensional array #' -#' @param df a [data.frame()] (or [tibble::tibble()]) in long format (each row contains one value for the intended array) -#' @param dim_var_names a [character()] vector of variable names in `df`. All of these variables should be factors, or a warning will be produced. -#' @param value_var_name a [character()] variable containing a variable name from `df` which contains the values for the intended array. -#' @return an [array()] with dimensions defined by the variables in `df` listed in `dim_var_names` +#' @param df a [data.frame()] (or [tibble::tibble()]) in long format +#' (each row contains one value for the intended array) +#' @param dim_var_names a [character()] vector of variable names in `df`. +#' All of these variables should be factors, or a warning will be produced. +#' @param value_var_name a [character()] variable containing a variable name +#' from `df` which contains the values for the intended array. +#' @return an [array()] with dimensions defined by the variables in `df` +#' listed in `dim_var_names` #' #' @examples #' library(dplyr) @@ -33,9 +37,11 @@ df.to.array <- function( #' cols = c("Sepal.Length", "Sepal.Width", "Petal.Width", "Petal.Length") #' ) %>% #' mutate(parameter = factor(parameter, levels = unique(parameter))) -#' arr <- df %>% serocalculator:::df.to.array(dim_var_names = c("parameter", "Species")) +#' arr <- df %>% +#' serocalculator:::df_to_array( +#' dim_var_names = c("parameter", "Species")) #' ftable(arr[,,1:5]) -#' @noRd +#' @dev df_to_array <- function( df, dim_var_names, @@ -49,7 +55,7 @@ df_to_array <- function( all_factors <- df %>% - select(dim_var_names) %>% + select(all_of(dim_var_names)) %>% sapply(FUN = is.factor) %>% all() @@ -77,6 +83,6 @@ df_to_array <- function( ) df %>% - mutate(.by = all_of(dim_var_names), obs = 1:n()) %>% + mutate(.by = all_of(dim_var_names), obs = row_number()) %>% xtabs(formula = formula(xtabs_formula)) } diff --git a/R/est.incidence.R b/R/est.incidence.R index fb9c8dc8..a0c291e5 100644 --- a/R/est.incidence.R +++ b/R/est.incidence.R @@ -119,7 +119,7 @@ est.incidence <- function( if (build_graph) { if (verbose) message("building likelihood graph") - graph <- graph.loglik( + graph <- graph_loglik( highlight_points = lambda_start, highlight_point_names = "lambda_start", pop_data = pop_data, diff --git a/R/graph.decay.curves.R b/R/graph.decay.curves.R index 963bde5d..a9248be8 100644 --- a/R/graph.decay.curves.R +++ b/R/graph.decay.curves.R @@ -4,22 +4,32 @@ #' @param verbose verbose output #' @param xlim range of x values to graph #' @param n_curves how many curves to plot (see details). -#' @param n_points Number of points to interpolate along the x axis (passed to [ggplot2::geom_function()]) -#' @param rows_to_graph which rows of `curve_params` to plot (overrides `n_curves`). -#' @param alpha (passed to [ggplot2::geom_function()]) how transparent the curves should be: +#' @param n_points Number of points to interpolate along the x axis +#' (passed to [ggplot2::geom_function()]) +#' @param rows_to_graph which rows of `curve_params` to plot +#' (overrides `n_curves`). +#' @param alpha (passed to [ggplot2::geom_function()]) +#' how transparent the curves should be: #' * 0 = fully transparent (invisible) #' * 1 = fully opaque -#' @param log_x should the x-axis be on a logarithmic scale (`TRUE`) or linear scale (`FALSE`, default)? -#' @param log_y should the Y-axis be on a logarithmic scale (default, `TRUE`) or linear scale (`FALSE`)? +#' @param log_x should the x-axis be on a logarithmic scale (`TRUE`) +#' or linear scale (`FALSE`, default)? +#' @param log_y should the Y-axis be on a logarithmic scale +#' (default, `TRUE`) or linear scale (`FALSE`)? #' @inheritParams ggplot2::geom_function #' @inheritDotParams ggplot2::geom_function #' @returns a [ggplot2::ggplot()] object #' @details #' ## `n_curves` and `rows_to_graph` -#' In most cases, `curve_params` will contain too many rows of MCMC samples for all of these samples to be plotted at once. -#' * Setting the `n_curves` argument to a value smaller than the number of rows in `curve_params` will cause this function to select the first `n_curves` rows to graph. -#' * Setting `n_curves` larger than the number of rows in ` will result all curves being plotted. -#' * If the user directly specifies the `rows_to_graph` argument, then `n_curves` has no effect. +#' In most cases, `curve_params` will contain too many rows of MCMC +#' samples for all of these samples to be plotted at once. +#' * Setting the `n_curves` argument to a value smaller than the +#' number of rows in `curve_params` will cause this function to select +#' the first `n_curves` rows to graph. +#' * Setting `n_curves` larger than the number of rows in ` will +#' result all curves being plotted. +#' * If the user directly specifies the `rows_to_graph` argument, +#' then `n_curves` has no effect. #' @examples #' library(dplyr) # loads the `%>%` operator and `dplyr::filter()` #' @@ -27,7 +37,7 @@ #' dplyr::filter(antigen_iso == "HlyE_IgG") %>% #' serocalculator:::plot_curve_params_one_ab() #' -plot_curve_params_one_ab = function( +plot_curve_params_one_ab <- function( object, verbose = FALSE, alpha = .4, @@ -35,44 +45,30 @@ plot_curve_params_one_ab = function( n_points = 1000, log_x = FALSE, log_y = TRUE, - rows_to_graph = 1:min(n_curves, nrow(object)), - xlim = c(10^-1, 10^3.1), + rows_to_graph = seq_len(min(n_curves, nrow(object))), + xlim = c(10 ^ -1, 10 ^ 3.1), ...) { plot1 <- ggplot2::ggplot() + # ggplot2::scale_x_log10() + ggplot2::theme_linedraw() + - ggplot2::theme( - axis.line = ggplot2::element_line() - ) + - ggplot2::labs( - x = "Days since fever onset", - y = "Antibody Concentration" - ) + - ggplot2::ggtitle("Decay Curve") + - ggplot2::theme( - plot.title = - ggplot2::element_text( - size = 20, - face = "bold" - ) - ) + ggplot2::theme(axis.line = ggplot2::element_line()) + + ggplot2::labs(x = "Days since fever onset", y = "Antibody concentration") + + ggplot2::ggtitle("Antibody Response Curve") + + ggplot2::theme(plot.title = + ggplot2::element_text(size = 20, face = "bold")) if (log_y) { - plot1 <- plot1 + - ggplot2::scale_y_log10( - # limits = c(0.9, 2000), - labels = scales::label_comma(), - # breaks = c(0.1, 1, 10, 100, 1000), - minor_breaks = NULL - ) + plot1 <- + plot1 + + ggplot2::scale_y_log10(labels = scales::label_comma(), + minor_breaks = NULL) } layer_function <- function(cur_row) { cur_params <- object[cur_row, ] ggplot2::geom_function( alpha = alpha, - # aes(color = cur_row), fun = ab0, args = list(curve_params = cur_params), n = n_points @@ -80,33 +76,16 @@ plot_curve_params_one_ab = function( } layers <- - lapply( - X = rows_to_graph, - FUN = layer_function - ) + lapply(X = rows_to_graph, FUN = layer_function) plot1 <- plot1 + layers if (log_x) { plot1 <- plot1 + - ggplot2::scale_x_log10( - limits = xlim, - labels = scales::label_comma() - ) + ggplot2::scale_x_log10(limits = xlim, labels = scales::label_comma()) } else { plot1 <- plot1 + ggplot2::xlim(xlim) } return(plot1) } - - -# ggplot() + -# geom_line(data = serocourse.all, aes(x= t, y = res, group = iter)) + -# facet_wrap(~antigen_iso, ncol=2) + -# scale_y_log10(limits = c(0.9, 2000), breaks = c(1, 10, 100, 1000), minor_breaks = NULL) + -# theme_minimal() + -# theme(axis.line=element_line()) + -# labs(x="Days since fever onset", y="ELISA units") - -# mcmc %>% ungroup() %>% slice_head(by = antigen_iso, n = 10) %>% droplevels() %>% plot_curve_params_one_ab(alpha = .4) %>% print() diff --git a/R/graph.loglik.R b/R/graph.loglik.R index a8d5b699..7e97d8d4 100644 --- a/R/graph.loglik.R +++ b/R/graph.loglik.R @@ -3,8 +3,10 @@ #' @param highlight_points a possible highlighted value #' @param x sequence of lambda values to graph #' @param highlight_point_names labels for highlighted points -#' @param log_x should the x-axis be on a logarithmic scale (`TRUE`) or linear scale (`FALSE`, default)? -#' @param previous_plot if not NULL, the current data is added to the existing graph +#' @param log_x should the x-axis be on a logarithmic scale (`TRUE`) +#' or linear scale (`FALSE`, default)? +#' @param previous_plot if not NULL, the current data is added to the +#' existing graph #' @param curve_label if not NULL, add a label for the curve #' @inheritParams log_likelihood #' @inheritDotParams log_likelihood -lambda @@ -33,21 +35,22 @@ #' y.high = c(5e6, 5e6)) # High cutoff (y.high) #' #' # Graph the log likelihood -#' lik_HlyE_IgA <- graph.loglik( -#' pop_data = xs_data, -#' curve_params = dmcmc, -#' noise_params = cond, -#' antigen_isos = "HlyE_IgA", -#' log_x = TRUE +#' lik_HlyE_IgA <- # nolint: object_name_linter +#' graph_loglik( +#' pop_data = xs_data, +#' curve_params = dmcmc, +#' noise_params = cond, +#' antigen_isos = "HlyE_IgA", +#' log_x = TRUE #' ) #' -#' lik_HlyE_IgA +#' lik_HlyE_IgA # nolint: object_name_linter #' -graph.loglik = function( +graph_loglik <- function( pop_data, curve_params, noise_params, - antigen_isos, + antigen_isos = pop_data %>% get_biomarker_levels(), x = 10^seq(-3, 0, by = .1), highlight_points = NULL, highlight_point_names = "highlight_points", @@ -55,8 +58,11 @@ graph.loglik = function( previous_plot = NULL, curve_label = paste(antigen_isos, collapse = " + "), ...) { - if (!is.list(curve_params) && - !is.element("d", names(curve_params))) { + + needs_rescale <- + !is.list(curve_params) && + !is.element("d", names(curve_params)) + if (needs_rescale) { curve_params <- curve_params %>% dplyr::mutate( diff --git a/R/load_pop_data.R b/R/load_pop_data.R index 8230df89..c1d711ea 100644 --- a/R/load_pop_data.R +++ b/R/load_pop_data.R @@ -1,6 +1,7 @@ #' Load a cross-sectional antibody survey data set #' -#' @param file_path path to an RDS file containing a cross-sectional antibody survey data set, stored as a [data.frame()] or [tibble::tbl_df] +#' @param file_path path to an RDS file containing a cross-sectional antibody +#' survey data set, stored as a [data.frame()] or [tibble::tbl_df] #' @inheritDotParams as_pop_data #' @returns a `pop_data` object (a [tibble::tbl_df] with extra attributes) #' @export @@ -8,9 +9,7 @@ #' xs_data <- load_pop_data("https://osf.io/download//n6cp3/") #' #' print(xs_data) -load_pop_data <- function(file_path, - ...) { - +load_pop_data <- function(file_path, ...) { pop_data <- file_path %>% readr::read_rds() %>% @@ -18,287 +17,3 @@ load_pop_data <- function(file_path, return(pop_data) } - -get_age <- function(object, ...) { - UseMethod("get_age", object) -} - -#' @export -get_age.pop_data <- function(object, ...) { - # get age data - age_data <- object %>% pull(attr(object, "age_var")) - - return(age_data) -} - -get_age_var <- function(object, ...) { - UseMethod("get_age_var", object) -} - -#' @export -get_age_var.pop_data <- function(object, ...) { - # get value attribute - age_var <- attributes(object)$age_var - - return(age_var) -} - -get_value <- function(object, ...) { - UseMethod("get_value", object) -} - -#' @export -get_value.pop_data <- function(object, ...) { - # get age data - value_data <- object %>% pull(attr(object, "value_var")) - - return(value_data) -} - -get_value_var <- function(object, ...) { - UseMethod("get_value_var", object) -} - -#' @export -get_value_var.pop_data <- function(object, ...) { - # get value attribute - value_var <- attributes(object)$value_var - - return(value_var) -} - -get_id <- function(object, ...) { - UseMethod("get_id", object) -} - -#' @export -get_id.pop_data <- function(object, ...) { - # get age data - id_data <- object %>% pull(attr(object, "id_var")) - - return(id_data) -} - -get_id_var <- function(object, ...) { - UseMethod("get_id_var", object) -} - -#' @export -get_id_var.pop_data <- function(object, ...) { - # get value attribute - id_var <- attributes(object)$id_var - - return(id_var) -} - -set_biomarker_var <- function(object, ...) { - UseMethod("set_biomarker_var", object) -} - -#' @export -set_biomarker_var.pop_data = function(object, - biomarker = "antigen_iso", - standardize = TRUE, - ...) -{ - if (biomarker %in% colnames(object)) - { - attr(object, "biomarker_var") <- biomarker - } else - { - cli::cli_abort('data does not include column "{biomarker}"') - } - - if (standardize) - { - object <- object %>% - rename(c("antigen_iso" = attr(object, "biomarker_var"))) - - # update attribute - attr(object, "biomarker_var") <- "antigen_iso" - } - - return(object) - -} - -get_biomarker_levels <- function(object, ...) -{ - UseMethod("get_biomarker_levels", object) -} - -#' @export -get_biomarker_levels.pop_data <- function(object, ...) -{ - attr(object, "antigen_isos") -} - -get_biomarker_names <- function(object, ...) { - UseMethod("get_biomarker_names", object) -} - -#' @export -get_biomarker_names.pop_data <- function(object, ...) { - # get biomarker name data - biomarker_data <- object %>% pull(get_biomarker_names_var(object)) - - return(biomarker_data) -} - -get_biomarker_names_var <- function(object, ...) { - UseMethod("get_biomarker_names_var", object) -} - -#' @export -get_biomarker_names_var.pop_data <- function(object, ...) { - # get value attribute - biomarker_var <- attributes(object)[["biomarker_var"]] - - return(biomarker_var) -} - - -set_age <- function(object, ...) { - UseMethod("set_age", object) -} - -#' @export -set_age.pop_data <- function(object, age = "Age", standardize = TRUE, ...) { - # check if age column exists - if (age %in% colnames(object)) { - attr(object, "age_var") <- age - } else { - cli::cli_warn('The specified `age` column "{age}" does not exist.') - - # search age variable from object - age_var <- - grep( - x = colnames(object), - value = TRUE, - pattern = age, - ignore.case = TRUE - ) - - if (length(age_var) == 1) { - attr(object, "age_var") <- age_var - - # create warning when using searched age instead of provided age - cli::cli_inform('Proceeding to use "{.var {age_var}}"') - } else if (length(age_var) == 0) { - cli::cli_abort("No similar column name was detected.") - } else if (length(age_var) > 1) - { - cli::cli_warn("Multiple potential matches found: {.var {age_var}}") - cli::cli_warn("Using first match: {.var {age_var[1]}}") - attr(object, "age_var") <- age_var[1] - } else - { - cli::cli_abort("{.code length(age_var)} = {.val {length(age_var)}}") - } - } - - if (standardize) { - object <- object %>% - rename(c("age" = attr(object, "age_var"))) - - # set age attribute - attr(object, "age_var") <- "age" - } - - return(object) -} - - -set_value <- function(object, ...) { - UseMethod("set_value", object) -} - -#' @export -set_value.pop_data <- function(object, value = "result", standardize = TRUE, ...) { - # check if value column exists - if (value %in% colnames(object)) { - attr(object, "value_var") <- value - } else { - cli::cli_warn('The specified `value` column "{.var {value}}" does not exist.') - - # search value variable from pop_data - value_var <- - grep( - x = colnames(object), - value = TRUE, - pattern = value, - ignore.case = TRUE - ) - - if (length(value_var) == 1) { - attr(object, "value_var") <- value_var - - # create warning when using searched age instead of provided age - cli::cli_inform('Proceeding to use "{.var {value_var}}"') - } else if (length(value_var) == 0) { - cli::cli_abort("No similar column name was detected.") - } else # if (length(value_var) > 1) - { - cli::cli_warn("Multiple potential matches found: {.var {value_var}}") - cli::cli_inform("Using first match: {.var {value_var[1]}}") - attr(object, "value_var") <- value_var[1] - } - } - - if (standardize) { - object <- object %>% - rename(c("value" = attr(object, "value_var"))) - - # set id attribute - attr(object, "value_var") <- "value" - } - - return(object) -} - -set_id <- function(object, ...) { - UseMethod("set_id", object) -} - -#' @export -set_id.pop_data <- function(object, id = "index_id", standardize = TRUE, ...) { - # check if id column exists - if (id %in% colnames(object)) { - attr(object, "id_var") <- id - } else { - cli::cli_warn('The specified {.var id} column {.val {id}} does not exist.') - - # search id variable from object - id_var <- - grep( - x = colnames(object), - value = TRUE, - pattern = id, - ignore.case = TRUE - ) - - if (length(id_var) == 1) { - attr(object, "id_var") <- id_var - - # create warning when using searched id instead of provided id - cli::cli_inform('Proceeding to use "{id_var}"') - } else if (length(id_var) == 0) { - cli::cli_abort("No similar column name was detected.") - } else # if (length(id_var) > 1) - { - cli::cli_warn("Multiple potential matches found: {.var {id_var}}") - cli::cli_inform("Using first match: {.var {id_var[1]}}") - attr(object, "id_var") <- id_var[1] - } - } - - if (standardize) { - object <- object %>% - rename(c("id" = attr(object, "id_var"))) - - # set id attribute - attr(object, "id_var") <- "id" - } - - return(object) -} diff --git a/R/llik.R b/R/log_likelihood.R similarity index 70% rename from R/llik.R rename to R/log_likelihood.R index 47b6c6d3..2022a83c 100644 --- a/R/llik.R +++ b/R/log_likelihood.R @@ -1,7 +1,9 @@ -#' Calculate log-likelihood +#' @title Calculate log-likelihood #' #' @description -#' Calculates the log-likelihood of a set of cross-sectional antibody response data, for a given incidence rate (`lambda`) value. +#' Calculates the log-likelihood of a set of cross-sectional antibody response +#' data, for a given incidence rate (`lambda`) value. +#' #' `r lifecycle::badge("deprecated")` #' #' `llik()` was renamed to [log_likelihood()] to create a more @@ -9,47 +11,45 @@ #' #' @keywords internal #' @export - -llik <- function( - lambda, - pop_data, - antigen_isos, - curve_params, - noise_params, - verbose = FALSE, - ...) { +llik <- function(...) { lifecycle::deprecate_warn("1.0.0", "llik()", "log_likelihood()") - log_likelihood(lambda, - pop_data, - antigen_isos, - curve_params, - noise_params, - verbose) + log_likelihood(...) } -#' Calculate log-likelihood +#' @title Calculate log-likelihood #' #' @description -#' Calculates the log-likelihood of a set of cross-sectional antibody response data, for a given incidence rate (`lambda`) value. -#' @param pop_data a [data.frame()] with cross-sectional serology data per antibody and age, and additional columns -#' @param antigen_isos Character vector listing one or more antigen isotypes. Values must match `pop_data`. -#' @param curve_params a [data.frame()] containing MCMC samples of parameters from the Bayesian posterior distribution of a longitudinal decay curve model. The parameter columns must be named: -#' - `antigen_iso`: a [character()] vector indicating antigen-isotype combinations +#' Calculates the log-likelihood of a set of cross-sectional antibody response +#' data, for a given incidence rate (`lambda`) value. +#' @param pop_data a [data.frame()] with cross-sectional serology data +#' by antibody and age, and additional columns +#' @param antigen_isos Character vector listing one or more antigen isotypes. +#' Values must match `pop_data`. +#' @param curve_params a [data.frame()] containing MCMC samples of parameters +#' from the Bayesian posterior distribution of a longitudinal decay curve model. +#' The parameter columns must be named: +#' - `antigen_iso`: a [character()] vector indicating antigen-isotype +#' combinations #' - `iter`: an [integer()] vector indicating MCMC sampling iterations #' - `y0`: baseline antibody level at $t=0$ ($y(t=0)$) #' - `y1`: antibody peak level (ELISA units) #' - `t1`: duration of infection -#' - `alpha`: antibody decay rate (1/days for the current longitudinal parameter sets) +#' - `alpha`: antibody decay rate +#' (1/days for the current longitudinal parameter sets) #' - `r`: shape factor of antibody decay -#' @param noise_params a [data.frame()] (or [tibble::tibble()]) containing the following variables, specifying noise parameters for each antigen isotype: -#' * `antigen_iso`: antigen isotype whose noise parameters are being specified on each row +#' @param noise_params a [data.frame()] (or [tibble::tibble()]) +#' containing the following variables, +#' specifying noise parameters for each antigen isotype: +#' * `antigen_iso`: antigen isotype whose noise parameters are being specified +#' on each row #' * `nu`: biological noise #' * `eps`: measurement noise #' * `y.low`: lower limit of detection for the current antigen isotype #' * `y.high`: upper limit of detection for the current antigen isotype #' #' @param verbose logical: if TRUE, print verbose log information to console -#' @param ... additional arguments passed to other functions (not currently used). +#' @param ... additional arguments passed to other functions +#' (not currently used). #' @inheritParams f_dev #' @return the log-likelihood of the data with the current parameter values #' @export @@ -84,19 +84,19 @@ llik <- function( log_likelihood <- function( lambda, pop_data, - antigen_isos, curve_params, noise_params, + antigen_isos = get_biomarker_levels(pop_data), verbose = FALSE, ...) { # Start with zero total - nllTotal <- 0 + nll_total <- 0 # Loop over antigen_isos - for (cur_antibody in antigen_isos) - { + for (cur_antibody in antigen_isos) { # the inputs can be lists, after `split(~antigen_ios)` - # this gives some speedups compared to running filter() every time .nll() is called + # this gives some speedups compared to running filter() every time .nll() + # is called if (!is.data.frame(pop_data)) { cur_data <- pop_data[[cur_antibody]] cur_curve_params <- curve_params[[cur_antibody]] @@ -124,7 +124,7 @@ log_likelihood <- function( } } - nllSingle <- + nll_single <- f_dev( lambda = lambda, csdata = cur_data, @@ -133,11 +133,12 @@ log_likelihood <- function( ) # if (!is.na(nllSingle)) # not meaningful for vectorized f_dev() - { - nllTotal <- nllTotal + nllSingle # DEM note: summing log likelihoods represents an independence assumption for multiple Antibodies, given time since seroconversion - } + nll_total <- nll_total + nll_single + # note: summing log likelihoods represents an independence assumption + # for multiple Antibodies, given time since seroconversion + } # Return total log-likelihood - return(-nllTotal) + return(-nll_total) } diff --git a/R/reexports.R b/R/reexports.R new file mode 100644 index 00000000..88758ee8 --- /dev/null +++ b/R/reexports.R @@ -0,0 +1,3 @@ +#' @export +#' @importFrom ggplot2 autoplot +ggplot2::autoplot diff --git a/R/sees_baseline_data.R b/R/sees_baseline_data.R deleted file mode 100644 index 90f878da..00000000 --- a/R/sees_baseline_data.R +++ /dev/null @@ -1,14 +0,0 @@ -#' #' @docType data -#' #' -#' #' @name sees_crossSectional_baseline_allCountries -#' #' -#' #' @title -#' #' Cross-sectional Typhoid data -#' #' -#' #' @description -#' #' A [tibble::tibble()] -#' #' -#' #' @usage -#' #' sees_crossSectional_baseline_allCountries -#' #' -#' NULL diff --git a/R/serocalculator-package.R b/R/serocalculator-package.R index 40da54fb..b1966566 100644 --- a/R/serocalculator-package.R +++ b/R/serocalculator-package.R @@ -8,29 +8,10 @@ #' Estimating Infection Rates from Serological Data #' #' @description -#' This package translates antibody levels measured in a (cross-sectional) population sample into an -#' estimate of the frequency with which seroconversions (infections) occur in the sampled population. -#' -#' The API for this package includes the following functions: -#' -#' Function Name | Purpose -#' ----------------- | -------------------------------------------- -#' [load_pop_data()] | loading cross-sectional antibody survey data -#' [check_pop_data()] | checking antibody data -#' [summary.pop_data()] | numerical summaries of antibody data -#' [autoplot.pop_data()] | graphs of antibody data distributions -#' [load_curve_params()] | loading antibody decay curve models -#' [autoplot.curve_params()] | graphing antibody decay curves -#' [llik()] | computing log-likelihoods -#' [graph.loglik()] | graphing log-likelihood functions -#' [autoplot.seroincidence()] | graphing log-likelihood functions -#' [autoplot.seroincidence.by()] | graphing log-likelihood functions -#' [est.incidence()] | estimating incidence rates -#' [est.incidence.by()] | estimating incidence rates by strata -#' [summary.seroincidence.by()] | summarizing stratified incidence rate estimates -#' [autoplot.summary.seroincidence.by()] | graphing incidence rate estimates -#' [sim.cs()] | simulating cross-sectional population antibody data using longitudinal seroresponse models -#' +#' This package translates antibody levels measured in a +#' (cross-sectional) population sample into an +#' estimate of the frequency with which seroconversions (infections) +#' occur in the sampled population. #' #' @author #' * Peter Teunis \email{p.teunis@@emory.edu} @@ -42,59 +23,99 @@ #' #' ***Methods for estimating seroincidence*** #' -#' * Teunis, P. F. M., and J. C. H. van Eijkeren. "Estimation of seroconversion rates for infectious diseases: Effects of age and noise." Statistics in Medicine 39, no. 21 (2020): 2799-2814. +#' * Teunis, P. F. M., and J. C. H. van Eijkeren. +#' "Estimation of seroconversion rates for infectious diseases: +#' Effects of age and noise." +#' Statistics in Medicine 39, no. 21 (2020): 2799-2814. #' -#' * Teunis, P. F. M., J. C. H. van Eijkeren, W. F. de Graaf, A. Bonačić Marinović, and M. E. E. Kretzschmar. "Linking the seroresponse to infection to within-host heterogeneity in antibody production." Epidemics 16 (2016): 33-39. +#' * Teunis, P. F. M., J. C. H. van Eijkeren, W. F. de Graaf, +#' A. Bonačić Marinović, and M. E. E. Kretzschmar. +#' "Linking the seroresponse to infection to within-host heterogeneity +#' in antibody production." Epidemics 16 (2016): 33-39. #' #' #' ***Applications*** #' -#' * Aiemjoy, Kristen, Jessica C. Seidman, Senjuti Saha, Sira Jam Munira, Mohammad Saiful Islam Sajib, Syed Muktadir Al Sium, Anik Sarkar et al. "Estimating typhoid incidence from community-based serosurveys: a multicohort study." The Lancet Microbe 3, no. 8 (2022): e578-e587. -#' -#' * Aiemjoy, Kristen, John Rumunu, Juma John Hassen, Kirsten E. Wiens, Denise Garrett, Polina Kamenskaya, Jason B. Harris et al. "Seroincidence of enteric fever, Juba, South Sudan." Emerging infectious diseases 28, no. 11 (2022): 2316. -#' -#' * Monge, S., Teunis, P. F., Friesema, I., Franz, E., Ang, W., van Pelt, W., Mughini-Gras, L. -#' "Immune response-eliciting exposure to Campylobacter vastly exceeds the incidence of clinically -#' overt campylobacteriosis but is associated with similar risk factors: A nationwide serosurvey in the Netherlands" +#' * Aiemjoy, Kristen, Jessica C. Seidman, Senjuti Saha, +#' Sira Jam Munira, Mohammad Saiful Islam Sajib, +#' Syed Muktadir Al Sium, Anik Sarkar et al. +#' "Estimating typhoid incidence from community-based serosurveys: +#' a multicohort study." The Lancet Microbe 3, no. 8 (2022): e578-e587. +#' +#' * Aiemjoy, Kristen, John Rumunu, Juma John Hassen, +#' Kirsten E. Wiens, Denise Garrett, Polina Kamenskaya, +#' Jason B. Harris et al. "Seroincidence of enteric fever, Juba, +#' South Sudan." Emerging infectious diseases 28, no. 11 (2022): 2316. +#' +#' * Monge, S., Teunis, P. F., Friesema, I., Franz, E., Ang, W., +#' van Pelt, W., Mughini-Gras, L. +#' "Immune response-eliciting exposure to Campylobacter vastly exceeds +#' the incidence of clinically +#' overt campylobacteriosis but is associated with similar +#' risk factors: A nationwide serosurvey in the Netherlands" #' Journal of Infection, 2018, 1--7, doi:10.1016/j.jinf.2018.04.016 #' #' * Kretzschmar, M., Teunis, P. F., Pebody, R. G. -#' "Incidence and reproduction numbers of pertussis: estimates from serological and social contact data in five European countries" -#' PLoS Medicine 7, no. 6 (June 1, 2010):e1000291. doi:10.1371/journal.pmed.1000291. -#' -#' * Simonsen, J., Strid, M. A., Molbak, K., Krogfelt, K. A., Linneberg, A., Teunis, P. -#' "Sero-epidemiology as a tool to study the incidence of Salmonella infections in humans" -#' Epidemiology and Infection 136, no. 7 (July 1, 2008): 895--902. doi:10.1017/S0950268807009314 -#' -#' * Simonsen, J., Teunis, P. F., van Pelt, W., van Duynhoven, Y., Krogfelt, K. A., Sadkowska-Todys, M., Molbak K. -#' "Usefulness of seroconversion rates for comparing infection pressures between countries" -#' Epidemiology and Infection, April 12, 2010, 1--8. doi:10.1017/S0950268810000750. -#' -#' * Falkenhorst, G., Simonsen, J., Ceper, T. H., van Pelt, W., de Valk, H., Sadkowska-Todys, M., Zota, L., Kuusi, M., Jernberg, C., Rota, M. C., van Duynhoven, Y. T., Teunis, P. F., Krogfelt, K. A., Molbak, K. -#' "Serological cross-sectional studies on salmonella incidence in eight European countries: no correlation with incidence of reported cases" -#' BMC Public Health 12, no. 1 (July 15, 2012): 523--23. doi:10.1186/1471-2458-12-523. -#' -#' * Teunis, P. F., Falkenhorst, G., Ang, C. W., Strid, M. A., De Valk, H., Sadkowska-Todys, M., Zota, L., Kuusi, M., Rota, M. C., Simonsen, J. B., Molbak, K., Van Duynhoven, Y. T., van Pelt, W. +#' "Incidence and reproduction numbers of pertussis: estimates from +#' serological and social contact data in five European countries" +#' PLoS Medicine 7, no. 6 (June 1, 2010):e1000291. +#' doi:10.1371/journal.pmed.1000291. +#' +#' * Simonsen, J., Strid, M. A., Molbak, K., Krogfelt, K. A., +#' Linneberg, A., Teunis, P. +#' "Sero-epidemiology as a tool to study the incidence of +#' Salmonella infections in humans" +#' Epidemiology and Infection 136, no. 7 (July 1, 2008): 895--902. +#' doi:10.1017/S0950268807009314 +#' +#' * Simonsen, J., Teunis, P. F., van Pelt, W., van Duynhoven, Y., +#' Krogfelt, K. A., Sadkowska-Todys, M., Molbak K. +#' "Usefulness of seroconversion rates for comparing infection +#' pressures between countries" +#' Epidemiology and Infection, April 12, 2010, 1--8. +#' doi:10.1017/S0950268810000750. +#' +#' * Falkenhorst, G., Simonsen, J., Ceper, T. H., van Pelt, W., +#' de Valk, H., Sadkowska-Todys, M., Zota, L., Kuusi, M., +#' Jernberg, C., Rota, M. C., van Duynhoven, Y. T., Teunis, P. F., +#' Krogfelt, K. A., Molbak, K. +#' "Serological cross-sectional studies on salmonella incidence in +#' eight European countries: no correlation with incidence of +#' reported cases" +#' BMC Public Health 12, no. 1 (July 15, 2012): 523--23. +#' doi:10.1186/1471-2458-12-523. +#' +#' * Teunis, P. F., Falkenhorst, G., Ang, C. W., Strid, M. A., +#' De Valk, H., Sadkowska-Todys, M., Zota, L., Kuusi, M., Rota, M. C., Simonsen, J. B., Molbak, K., Van Duynhoven, Y. T., van Pelt, W. #' "Campylobacter seroconversion rates in selected countries in the European Union" -#' Epidemiology and Infection 141, no. 10 (2013): 2051--57. doi:10.1017/S0950268812002774. +#' Epidemiology and Infection 141, no. 10 (2013): 2051--57. +#' doi:10.1017/S0950268812002774. #' -#' * de Melker, H. E., Versteegh, F. G., Schellekens, J. F., Teunis, P. F., Kretzschmar, M. -#' "The incidence of Bordetella pertussis infections estimated in the population from a combination of serological surveys" +#' * de Melker, H. E., Versteegh, F. G., Schellekens, J. F., +#' Teunis, P. F., Kretzschmar, M. +#' "The incidence of Bordetella pertussis infections estimated in the +#' population from a combination of serological surveys" #' The Journal of Infection 53, no. 2 (August 1, 2006): 106--13. doi:10.1016/j.jinf.2005.10.020 #' #' #' ***Quantification of seroresponse*** #' #' * de Graaf, W. F., Kretzschmar, M. E., Teunis, P. F., Diekmann, O. -#' "A two-phase within-host model for immune response and its application to serological profiles of pertussis" +#' "A two-phase within-host model for immune response and its +#' application to serological profiles of pertussis" #' Epidemics 9 (2014):1--7. doi:10.1016/j.epidem.2014.08.002. #' -#' * Berbers, G. A., van de Wetering, M. S., van Gageldonk, P. G., Schellekens, J. F., Versteegh, F. G., Teunis, P.F. -#' "A novel method for evaluating natural and vaccine induced serological responses to Bordetella pertussis antigens" -#' Vaccine 31, no. 36 (August 12, 2013): 3732--38. doi:10.1016/j.vaccine.2013.05.073. -#' -#' * Versteegh, F. G., Mertens, P. L., de Melker, H. E., Roord, J. J., Schellekens, J. F., Teunis, P. F. -#' "Age-specific long-term course of IgG antibodies to pertussis toxin after symptomatic infection with Bordetella pertussis" +#' * Berbers, G. A., van de Wetering, M. S., van Gageldonk, P. G., +#' Schellekens, J. F., Versteegh, F. G., Teunis, P.F. +#' "A novel method for evaluating natural and vaccine induced +#' serological responses to Bordetella pertussis antigens" +#' Vaccine 31, no. 36 (August 12, 2013): 3732--38. +#' doi:10.1016/j.vaccine.2013.05.073. +#' +#' * Versteegh, F. G., Mertens, P. L., de Melker, H. E., Roord, J. J., +#' Schellekens, J. F., Teunis, P. F. +#' "Age-specific long-term course of IgG antibodies to pertussis toxin +#' after symptomatic infection with Bordetella pertussis" #' Epidemiology and Infection 133, no. 4 (August 1, 2005): 737--48. #' #' * Teunis, P. F., van der Heijden, O. G., de Melker, H. E., Schellekens, J. F., Versteegh, F. G., Kretzschmar, M. E. diff --git a/R/set_biomarker_var.R b/R/set_biomarker_var.R new file mode 100644 index 00000000..3aead4d2 --- /dev/null +++ b/R/set_biomarker_var.R @@ -0,0 +1,22 @@ +set_biomarker_var <- function(object, + biomarker = "antigen_iso", + standardize = TRUE, + ...) { + if (biomarker %in% colnames(object)) { + attr(object, "biomarker_var") <- biomarker + } else { + cli::cli_abort('data does not include column "{biomarker}"', + class = "missing variable") + } + + if (standardize) { + object <- object %>% + rename(c("antigen_iso" = attr(object, "biomarker_var"))) + + # update attribute + attr(object, "biomarker_var") <- "antigen_iso" + } + + return(object) + +} diff --git a/README.md b/README.md index 99d73bd6..48e4fe84 100644 --- a/README.md +++ b/README.md @@ -11,6 +11,10 @@ [![codecov](https://codecov.io/gh/UCD-SERG/serocalculator/graph/badge.svg?token=85CXV6GN2T)](https://codecov.io/gh/UCD-SERG/serocalculator) + #> Warning: package 'qrcode' was built under R version 4.4.1 + +![QR code for {serocalculator} website](man/figures/qr.svg) + Antibody levels measured in a cross–sectional population sample can be translated into an estimate of the frequency with which seroconversions (infections) occur in the sampled population. In other words, the diff --git a/README.qmd b/README.qmd index 48162420..179427b5 100644 --- a/README.qmd +++ b/README.qmd @@ -1,7 +1,6 @@ --- title: "serocalculator" format: gfm -bibliography: vignettes/articles/references.bib --- @@ -22,6 +21,15 @@ knitr::opts_chunk$set( [![codecov](https://codecov.io/gh/UCD-SERG/serocalculator/graph/badge.svg?token=85CXV6GN2T)](https://codecov.io/gh/UCD-SERG/serocalculator) +```{r svg} +#| echo: false +#| label: fig-qr-code +library(qrcode) +code <- qr_code("https://ucd-serg.github.io/serocalculator/") +generate_svg(code, filename = "man/figures/qr.svg") +``` + +![QR code for {serocalculator} website](man/figures/qr.svg) Antibody levels measured in a cross–sectional population sample can be translated into an estimate of the frequency with which seroconversions (infections) occur in the sampled population. In other words, the presence of many high antibody titers indicates that many individuals likely experienced infection recently and the burden of disease is high in the population, diff --git a/_quarto.yml b/_quarto.yml new file mode 100644 index 00000000..cd0da652 --- /dev/null +++ b/_quarto.yml @@ -0,0 +1,29 @@ +project: + render: + - '*.qmd' +author: "UC Davis Seroepidemiology Research Group (UCD-SERG)" +date: '`r Sys.Date()`' +format: + html: + number-sections: true + toc: true + html-math-method: mathjax + fig-cap-location: top + echo: false + revealjs: + reference-location: document + echo: false + html-math-method: mathjax + scrollable: true + auto-stretch: true + fig-height: 3 + fig-width: 7 + number-sections: true + number-depth: 2 + fig-cap-location: top + pptx: + number-sections: false + slide-level: 2 + reference-location: document + # reference-doc: vignettes/articles/templates/Gates_template.pptx +bibliography: vignettes/references.bib diff --git a/data-raw/sees_baseline_data.R b/data-raw/sees_baseline_data.R new file mode 100644 index 00000000..472c0e6e --- /dev/null +++ b/data-raw/sees_baseline_data.R @@ -0,0 +1,14 @@ +# #' @docType data +# #' +# #' @name sees_crossSectional_baseline_allCountries +# #' +# #' @title +# #' Cross-sectional Typhoid data +# #' +# #' @description +# #' A [tibble::tibble()] +# #' +# #' @usage +# #' sees_crossSectional_baseline_allCountries +# #' +# NULL diff --git a/inst/WORDLIST b/inst/WORDLIST index 49fbc668..f03a07b0 100644 --- a/inst/WORDLIST +++ b/inst/WORDLIST @@ -1,5 +1,3 @@ -95th -AKU Aga Ang Anik @@ -20,7 +18,6 @@ Falkenhorst Fogarty Friesema Gageldonk -GitHub Graaf Gras Grolemund @@ -34,12 +31,14 @@ IgM Jernberg Juba Juma -KGH Kamenskaya Kretzschmar Krogfelt Kuusi +LV +Lik Linneberg +Lotka Marinović Melker Mertens @@ -56,13 +55,11 @@ Orientia PLoS Pebody Polina -PRs RStudio Rai Roord Rtools Rumunu -SErologic Sadkowska Saha Saiful @@ -72,6 +69,7 @@ Schellekens Seidman Senjuti Sero +SeroEpidemiology Seroincidence Serological Simonsen @@ -81,48 +79,95 @@ Strid Syed TW Todys +Unif Valk -Vellore Vectorized +Vellore Versteegh +VignetteEncoding +VignetteEngine +VignetteIndexEntry +Volterra Wetering Wiens Zota al +ast behaviour +bioassays biomarker +boldsymbol +bookdown +callout campylobacteriosis +cdot +codecov de der devtools +displaystyle doi +dt epidem et +expf +forall +frac fx +geoms +ggproto hemolysin +infty +invf isos isotype isotypes jinf jitter kDa +knitr +le +leq +llik +logf +mathbb +mathbf +mathcal mc mcmc modelled multicohort +nd nr +overcount +overcounts +overline param params +pipetting pmed +pmf +pre +predation +providecommand +qmd +qquad recombinant +renewcommand +rescale +rmarkdown savePath +sectionally +sera seroconversion seroconversions +seroconverted seroepidemiology seroincidence seroincidences serologic serological +seronegative seroresponse seroresponses serosurvey @@ -133,11 +178,11 @@ subfigures th tibble titers +toc tsutsugamushi +undercount unstratified varepsilon -SeroEpidemiology -callout -codecov -geoms -ggproto +vec +vee +xsectionalData diff --git a/man/as_curve_params.Rd b/man/as_curve_params.Rd index 3f67702c..5bdb2cd7 100644 --- a/man/as_curve_params.Rd +++ b/man/as_curve_params.Rd @@ -9,10 +9,12 @@ as_curve_params(data, antigen_isos = NULL) \arguments{ \item{data}{a \code{\link[=data.frame]{data.frame()}} or \link[tibble:tbl_df-class]{tibble::tbl_df}} -\item{antigen_isos}{\code{\link[=character]{character()}} vector of antigen isotypes to be used in analyses} +\item{antigen_isos}{a \code{\link[=character]{character()}} vector of antigen isotypes +to be used in analyses} } \value{ -a \code{curve_data} object (a \link[tibble:tbl_df-class]{tibble::tbl_df} with extra attribute \code{antigen_isos}) +a \code{curve_data} object +(a \link[tibble:tbl_df-class]{tibble::tbl_df} with extra attribute \code{antigen_isos}) } \description{ Load antibody decay curve parameter diff --git a/man/autoplot.curve_params.Rd b/man/autoplot.curve_params.Rd index e32c4fde..933c45f6 100644 --- a/man/autoplot.curve_params.Rd +++ b/man/autoplot.curve_params.Rd @@ -24,15 +24,20 @@ \item{\code{verbose}}{verbose output} \item{\code{xlim}}{range of x values to graph} \item{\code{n_curves}}{how many curves to plot (see details).} - \item{\code{n_points}}{Number of points to interpolate along the x axis (passed to \code{\link[ggplot2:geom_function]{ggplot2::geom_function()}})} - \item{\code{rows_to_graph}}{which rows of \code{curve_params} to plot (overrides \code{n_curves}).} - \item{\code{alpha}}{(passed to \code{\link[ggplot2:geom_function]{ggplot2::geom_function()}}) how transparent the curves should be: + \item{\code{n_points}}{Number of points to interpolate along the x axis +(passed to \code{\link[ggplot2:geom_function]{ggplot2::geom_function()}})} + \item{\code{rows_to_graph}}{which rows of \code{curve_params} to plot +(overrides \code{n_curves}).} + \item{\code{alpha}}{(passed to \code{\link[ggplot2:geom_function]{ggplot2::geom_function()}}) +how transparent the curves should be: \itemize{ \item 0 = fully transparent (invisible) \item 1 = fully opaque }} - \item{\code{log_x}}{should the x-axis be on a logarithmic scale (\code{TRUE}) or linear scale (\code{FALSE}, default)?} - \item{\code{log_y}}{should the Y-axis be on a logarithmic scale (default, \code{TRUE}) or linear scale (\code{FALSE})?} + \item{\code{log_x}}{should the x-axis be on a logarithmic scale (\code{TRUE}) +or linear scale (\code{FALSE}, default)?} + \item{\code{log_y}}{should the Y-axis be on a logarithmic scale +(default, \code{TRUE}) or linear scale (\code{FALSE})?} }} } \value{ diff --git a/man/df_to_array.Rd b/man/df_to_array.Rd new file mode 100644 index 00000000..c915245f --- /dev/null +++ b/man/df_to_array.Rd @@ -0,0 +1,41 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/df_to_array.R +\name{df_to_array} +\alias{df_to_array} +\title{Convert a data.frame (or tibble) into a multidimensional array} +\usage{ +df_to_array(df, dim_var_names, value_var_name = "value") +} +\arguments{ +\item{df}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) in long format +(each row contains one value for the intended array)} + +\item{dim_var_names}{a \code{\link[=character]{character()}} vector of variable names in \code{df}. +All of these variables should be factors, or a warning will be produced.} + +\item{value_var_name}{a \code{\link[=character]{character()}} variable containing a variable name +from \code{df} which contains the values for the intended array.} +} +\value{ +an \code{\link[=array]{array()}} with dimensions defined by the variables in \code{df} +listed in \code{dim_var_names} +} +\description{ +Convert a data.frame (or tibble) into a multidimensional array +} +\examples{ +library(dplyr) +library(tidyr) + +df <- iris \%>\% + tidyr::pivot_longer( + names_to = "parameter", + cols = c("Sepal.Length", "Sepal.Width", "Petal.Width", "Petal.Length") + ) \%>\% + mutate(parameter = factor(parameter, levels = unique(parameter))) +arr <- df \%>\% + serocalculator:::df_to_array( + dim_var_names = c("parameter", "Species")) +ftable(arr[,,1:5]) +} +\keyword{internal} diff --git a/man/dot-nll.Rd b/man/dot-nll.Rd index 08ab530b..ba119ecf 100644 --- a/man/dot-nll.Rd +++ b/man/dot-nll.Rd @@ -12,21 +12,30 @@ \item{...}{ Arguments passed on to \code{\link[=log_likelihood]{log_likelihood}} \describe{ - \item{\code{pop_data}}{a \code{\link[=data.frame]{data.frame()}} with cross-sectional serology data per antibody and age, and additional columns} - \item{\code{antigen_isos}}{Character vector listing one or more antigen isotypes. Values must match \code{pop_data}.} - \item{\code{curve_params}}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters from the Bayesian posterior distribution of a longitudinal decay curve model. The parameter columns must be named: + \item{\code{pop_data}}{a \code{\link[=data.frame]{data.frame()}} with cross-sectional serology data +by antibody and age, and additional columns} + \item{\code{antigen_isos}}{Character vector listing one or more antigen isotypes. +Values must match \code{pop_data}.} + \item{\code{curve_params}}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters +from the Bayesian posterior distribution of a longitudinal decay curve model. +The parameter columns must be named: \itemize{ -\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype combinations +\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype +combinations \item \code{iter}: an \code{\link[=integer]{integer()}} vector indicating MCMC sampling iterations \item \code{y0}: baseline antibody level at $t=0$ ($y(t=0)$) \item \code{y1}: antibody peak level (ELISA units) \item \code{t1}: duration of infection -\item \code{alpha}: antibody decay rate (1/days for the current longitudinal parameter sets) +\item \code{alpha}: antibody decay rate +(1/days for the current longitudinal parameter sets) \item \code{r}: shape factor of antibody decay }} - \item{\code{noise_params}}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) containing the following variables, specifying noise parameters for each antigen isotype: + \item{\code{noise_params}}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) +containing the following variables, +specifying noise parameters for each antigen isotype: \itemize{ -\item \code{antigen_iso}: antigen isotype whose noise parameters are being specified on each row +\item \code{antigen_iso}: antigen isotype whose noise parameters are being specified +on each row \item \code{nu}: biological noise \item \code{eps}: measurement noise \item \code{y.low}: lower limit of detection for the current antigen isotype diff --git a/man/est.incidence.Rd b/man/est.incidence.Rd index 58b19797..ce207ed6 100644 --- a/man/est.incidence.Rd +++ b/man/est.incidence.Rd @@ -21,20 +21,27 @@ est.incidence( \arguments{ \item{pop_data}{a \link{data.frame} with cross-sectional serology data per antibody and age, and additional columns} -\item{curve_params}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters from the Bayesian posterior distribution of a longitudinal decay curve model. The parameter columns must be named: +\item{curve_params}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters +from the Bayesian posterior distribution of a longitudinal decay curve model. +The parameter columns must be named: \itemize{ -\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype combinations +\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype +combinations \item \code{iter}: an \code{\link[=integer]{integer()}} vector indicating MCMC sampling iterations \item \code{y0}: baseline antibody level at $t=0$ ($y(t=0)$) \item \code{y1}: antibody peak level (ELISA units) \item \code{t1}: duration of infection -\item \code{alpha}: antibody decay rate (1/days for the current longitudinal parameter sets) +\item \code{alpha}: antibody decay rate +(1/days for the current longitudinal parameter sets) \item \code{r}: shape factor of antibody decay }} -\item{noise_params}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) containing the following variables, specifying noise parameters for each antigen isotype: +\item{noise_params}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) +containing the following variables, +specifying noise parameters for each antigen isotype: \itemize{ -\item \code{antigen_iso}: antigen isotype whose noise parameters are being specified on each row +\item \code{antigen_iso}: antigen isotype whose noise parameters are being specified +on each row \item \code{nu}: biological noise \item \code{eps}: measurement noise \item \code{y.low}: lower limit of detection for the current antigen isotype diff --git a/man/est.incidence.by.Rd b/man/est.incidence.by.Rd index 123a4b16..3111521c 100644 --- a/man/est.incidence.by.Rd +++ b/man/est.incidence.by.Rd @@ -23,20 +23,27 @@ est.incidence.by( \arguments{ \item{pop_data}{a \link{data.frame} with cross-sectional serology data per antibody and age, and additional columns corresponding to each element of the \code{strata} input} -\item{curve_params}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters from the Bayesian posterior distribution of a longitudinal decay curve model. The parameter columns must be named: +\item{curve_params}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters +from the Bayesian posterior distribution of a longitudinal decay curve model. +The parameter columns must be named: \itemize{ -\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype combinations +\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype +combinations \item \code{iter}: an \code{\link[=integer]{integer()}} vector indicating MCMC sampling iterations \item \code{y0}: baseline antibody level at $t=0$ ($y(t=0)$) \item \code{y1}: antibody peak level (ELISA units) \item \code{t1}: duration of infection -\item \code{alpha}: antibody decay rate (1/days for the current longitudinal parameter sets) +\item \code{alpha}: antibody decay rate +(1/days for the current longitudinal parameter sets) \item \code{r}: shape factor of antibody decay }} -\item{noise_params}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) containing the following variables, specifying noise parameters for each antigen isotype: +\item{noise_params}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) +containing the following variables, +specifying noise parameters for each antigen isotype: \itemize{ -\item \code{antigen_iso}: antigen isotype whose noise parameters are being specified on each row +\item \code{antigen_iso}: antigen isotype whose noise parameters are being specified +on each row \item \code{nu}: biological noise \item \code{eps}: measurement noise \item \code{y.low}: lower limit of detection for the current antigen isotype diff --git a/man/figures/qr.svg b/man/figures/qr.svg new file mode 100644 index 00000000..40d5c446 --- /dev/null +++ b/man/figures/qr.svg @@ -0,0 +1,440 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/man/figures/qr_akm2.svg b/man/figures/qr_akm2.svg new file mode 100644 index 00000000..2e6e89cc --- /dev/null +++ b/man/figures/qr_akm2.svg @@ -0,0 +1,535 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/man/graph.loglik.Rd b/man/graph_loglik.Rd similarity index 73% rename from man/graph.loglik.Rd rename to man/graph_loglik.Rd index f4b51a69..efd15a33 100644 --- a/man/graph.loglik.Rd +++ b/man/graph_loglik.Rd @@ -1,14 +1,14 @@ % Generated by roxygen2: do not edit by hand % Please edit documentation in R/graph.loglik.R -\name{graph.loglik} -\alias{graph.loglik} +\name{graph_loglik} +\alias{graph_loglik} \title{Graph log-likelihood of data} \usage{ -graph.loglik( +graph_loglik( pop_data, curve_params, noise_params, - antigen_isos, + antigen_isos = pop_data \%>\% get_biomarker_levels(), x = 10^seq(-3, 0, by = 0.1), highlight_points = NULL, highlight_point_names = "highlight_points", @@ -19,29 +19,38 @@ graph.loglik( ) } \arguments{ -\item{pop_data}{a \code{\link[=data.frame]{data.frame()}} with cross-sectional serology data per antibody and age, and additional columns} +\item{pop_data}{a \code{\link[=data.frame]{data.frame()}} with cross-sectional serology data +by antibody and age, and additional columns} -\item{curve_params}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters from the Bayesian posterior distribution of a longitudinal decay curve model. The parameter columns must be named: +\item{curve_params}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters +from the Bayesian posterior distribution of a longitudinal decay curve model. +The parameter columns must be named: \itemize{ -\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype combinations +\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype +combinations \item \code{iter}: an \code{\link[=integer]{integer()}} vector indicating MCMC sampling iterations \item \code{y0}: baseline antibody level at $t=0$ ($y(t=0)$) \item \code{y1}: antibody peak level (ELISA units) \item \code{t1}: duration of infection -\item \code{alpha}: antibody decay rate (1/days for the current longitudinal parameter sets) +\item \code{alpha}: antibody decay rate +(1/days for the current longitudinal parameter sets) \item \code{r}: shape factor of antibody decay }} -\item{noise_params}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) containing the following variables, specifying noise parameters for each antigen isotype: +\item{noise_params}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) +containing the following variables, +specifying noise parameters for each antigen isotype: \itemize{ -\item \code{antigen_iso}: antigen isotype whose noise parameters are being specified on each row +\item \code{antigen_iso}: antigen isotype whose noise parameters are being specified +on each row \item \code{nu}: biological noise \item \code{eps}: measurement noise \item \code{y.low}: lower limit of detection for the current antigen isotype \item \code{y.high}: upper limit of detection for the current antigen isotype }} -\item{antigen_isos}{Character vector listing one or more antigen isotypes. Values must match \code{pop_data}.} +\item{antigen_isos}{Character vector listing one or more antigen isotypes. +Values must match \code{pop_data}.} \item{x}{sequence of lambda values to graph} @@ -49,9 +58,11 @@ graph.loglik( \item{highlight_point_names}{labels for highlighted points} -\item{log_x}{should the x-axis be on a logarithmic scale (\code{TRUE}) or linear scale (\code{FALSE}, default)?} +\item{log_x}{should the x-axis be on a logarithmic scale (\code{TRUE}) +or linear scale (\code{FALSE}, default)?} -\item{previous_plot}{if not NULL, the current data is added to the existing graph} +\item{previous_plot}{if not NULL, the current data is added to the +existing graph} \item{curve_label}{if not NULL, add a label for the curve} @@ -90,14 +101,15 @@ cond <- tibble( y.high = c(5e6, 5e6)) # High cutoff (y.high) # Graph the log likelihood -lik_HlyE_IgA <- graph.loglik( - pop_data = xs_data, - curve_params = dmcmc, - noise_params = cond, - antigen_isos = "HlyE_IgA", - log_x = TRUE +lik_HlyE_IgA <- # nolint: object_name_linter + graph_loglik( + pop_data = xs_data, + curve_params = dmcmc, + noise_params = cond, + antigen_isos = "HlyE_IgA", + log_x = TRUE ) -lik_HlyE_IgA +lik_HlyE_IgA # nolint: object_name_linter } diff --git a/man/llik.Rd b/man/llik.Rd index 38fa4961..4cf84bc8 100644 --- a/man/llik.Rd +++ b/man/llik.Rd @@ -1,21 +1,15 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/llik.R +% Please edit documentation in R/log_likelihood.R \name{llik} \alias{llik} \title{Calculate log-likelihood} \usage{ -llik( - lambda, - pop_data, - antigen_isos, - curve_params, - noise_params, - verbose = FALSE, - ... -) +llik(...) } \description{ -Calculates the log-likelihood of a set of cross-sectional antibody response data, for a given incidence rate (\code{lambda}) value. +Calculates the log-likelihood of a set of cross-sectional antibody response +data, for a given incidence rate (\code{lambda}) value. + \ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#deprecated}{\figure{lifecycle-deprecated.svg}{options: alt='[Deprecated]'}}}{\strong{[Deprecated]}} \code{llik()} was renamed to \code{\link[=log_likelihood]{log_likelihood()}} to create a more diff --git a/man/load_pop_data.Rd b/man/load_pop_data.Rd index 4dacd466..a2e87f95 100644 --- a/man/load_pop_data.Rd +++ b/man/load_pop_data.Rd @@ -7,7 +7,8 @@ load_pop_data(file_path, ...) } \arguments{ -\item{file_path}{path to an RDS file containing a cross-sectional antibody survey data set, stored as a \code{\link[=data.frame]{data.frame()}} or \link[tibble:tbl_df-class]{tibble::tbl_df}} +\item{file_path}{path to an RDS file containing a cross-sectional antibody +survey data set, stored as a \code{\link[=data.frame]{data.frame()}} or \link[tibble:tbl_df-class]{tibble::tbl_df}} \item{...}{ Arguments passed on to \code{\link[=as_pop_data]{as_pop_data}} diff --git a/man/log_likelihood.Rd b/man/log_likelihood.Rd index ade719ab..7a9267e9 100644 --- a/man/log_likelihood.Rd +++ b/man/log_likelihood.Rd @@ -1,5 +1,5 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/llik.R +% Please edit documentation in R/log_likelihood.R \name{log_likelihood} \alias{log_likelihood} \title{Calculate log-likelihood} @@ -7,9 +7,9 @@ log_likelihood( lambda, pop_data, - antigen_isos, curve_params, noise_params, + antigen_isos = get_biomarker_levels(pop_data), verbose = FALSE, ... ) @@ -17,39 +17,50 @@ log_likelihood( \arguments{ \item{lambda}{a \link{numeric} vector of incidence parameters, in events per person-year} -\item{pop_data}{a \code{\link[=data.frame]{data.frame()}} with cross-sectional serology data per antibody and age, and additional columns} +\item{pop_data}{a \code{\link[=data.frame]{data.frame()}} with cross-sectional serology data +by antibody and age, and additional columns} -\item{antigen_isos}{Character vector listing one or more antigen isotypes. Values must match \code{pop_data}.} - -\item{curve_params}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters from the Bayesian posterior distribution of a longitudinal decay curve model. The parameter columns must be named: +\item{curve_params}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters +from the Bayesian posterior distribution of a longitudinal decay curve model. +The parameter columns must be named: \itemize{ -\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype combinations +\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype +combinations \item \code{iter}: an \code{\link[=integer]{integer()}} vector indicating MCMC sampling iterations \item \code{y0}: baseline antibody level at $t=0$ ($y(t=0)$) \item \code{y1}: antibody peak level (ELISA units) \item \code{t1}: duration of infection -\item \code{alpha}: antibody decay rate (1/days for the current longitudinal parameter sets) +\item \code{alpha}: antibody decay rate +(1/days for the current longitudinal parameter sets) \item \code{r}: shape factor of antibody decay }} -\item{noise_params}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) containing the following variables, specifying noise parameters for each antigen isotype: +\item{noise_params}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) +containing the following variables, +specifying noise parameters for each antigen isotype: \itemize{ -\item \code{antigen_iso}: antigen isotype whose noise parameters are being specified on each row +\item \code{antigen_iso}: antigen isotype whose noise parameters are being specified +on each row \item \code{nu}: biological noise \item \code{eps}: measurement noise \item \code{y.low}: lower limit of detection for the current antigen isotype \item \code{y.high}: upper limit of detection for the current antigen isotype }} +\item{antigen_isos}{Character vector listing one or more antigen isotypes. +Values must match \code{pop_data}.} + \item{verbose}{logical: if TRUE, print verbose log information to console} -\item{...}{additional arguments passed to other functions (not currently used).} +\item{...}{additional arguments passed to other functions +(not currently used).} } \value{ the log-likelihood of the data with the current parameter values } \description{ -Calculates the log-likelihood of a set of cross-sectional antibody response data, for a given incidence rate (\code{lambda}) value. +Calculates the log-likelihood of a set of cross-sectional antibody response +data, for a given incidence rate (\code{lambda}) value. } \examples{ library(dplyr) diff --git a/man/plot_curve_params_one_ab.Rd b/man/plot_curve_params_one_ab.Rd index e5739f2d..9b312a07 100644 --- a/man/plot_curve_params_one_ab.Rd +++ b/man/plot_curve_params_one_ab.Rd @@ -12,7 +12,7 @@ plot_curve_params_one_ab( n_points = 1000, log_x = FALSE, log_y = TRUE, - rows_to_graph = 1:min(n_curves, nrow(object)), + rows_to_graph = seq_len(min(n_curves, nrow(object))), xlim = c(10^-1, 10^3.1), ... ) @@ -22,7 +22,8 @@ plot_curve_params_one_ab( \item{verbose}{verbose output} -\item{alpha}{(passed to \code{\link[ggplot2:geom_function]{ggplot2::geom_function()}}) how transparent the curves should be: +\item{alpha}{(passed to \code{\link[ggplot2:geom_function]{ggplot2::geom_function()}}) +how transparent the curves should be: \itemize{ \item 0 = fully transparent (invisible) \item 1 = fully opaque @@ -30,13 +31,17 @@ plot_curve_params_one_ab( \item{n_curves}{how many curves to plot (see details).} -\item{n_points}{Number of points to interpolate along the x axis (passed to \code{\link[ggplot2:geom_function]{ggplot2::geom_function()}})} +\item{n_points}{Number of points to interpolate along the x axis +(passed to \code{\link[ggplot2:geom_function]{ggplot2::geom_function()}})} -\item{log_x}{should the x-axis be on a logarithmic scale (\code{TRUE}) or linear scale (\code{FALSE}, default)?} +\item{log_x}{should the x-axis be on a logarithmic scale (\code{TRUE}) +or linear scale (\code{FALSE}, default)?} -\item{log_y}{should the Y-axis be on a logarithmic scale (default, \code{TRUE}) or linear scale (\code{FALSE})?} +\item{log_y}{should the Y-axis be on a logarithmic scale +(default, \code{TRUE}) or linear scale (\code{FALSE})?} -\item{rows_to_graph}{which rows of \code{curve_params} to plot (overrides \code{n_curves}).} +\item{rows_to_graph}{which rows of \code{curve_params} to plot +(overrides \code{n_curves}).} \item{xlim}{range of x values to graph} @@ -94,11 +99,16 @@ Graph an antibody decay curve model \details{ \subsection{\code{n_curves} and \code{rows_to_graph}}{ -In most cases, \code{curve_params} will contain too many rows of MCMC samples for all of these samples to be plotted at once. +In most cases, \code{curve_params} will contain too many rows of MCMC +samples for all of these samples to be plotted at once. \itemize{ -\item Setting the \code{n_curves} argument to a value smaller than the number of rows in \code{curve_params} will cause this function to select the first \code{n_curves} rows to graph. -\item Setting \code{n_curves} larger than the number of rows in ` will result all curves being plotted. -\item If the user directly specifies the \code{rows_to_graph} argument, then \code{n_curves} has no effect. +\item Setting the \code{n_curves} argument to a value smaller than the +number of rows in \code{curve_params} will cause this function to select +the first \code{n_curves} rows to graph. +\item Setting \code{n_curves} larger than the number of rows in ` will +result all curves being plotted. +\item If the user directly specifies the \code{rows_to_graph} argument, +then \code{n_curves} has no effect. } } } diff --git a/man/reexports.Rd b/man/reexports.Rd new file mode 100644 index 00000000..8570ff02 --- /dev/null +++ b/man/reexports.Rd @@ -0,0 +1,16 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/reexports.R +\docType{import} +\name{reexports} +\alias{reexports} +\alias{autoplot} +\title{Objects exported from other packages} +\keyword{internal} +\description{ +These objects are imported from other packages. Follow the links +below to see their documentation. + +\describe{ + \item{ggplot2}{\code{\link[ggplot2]{autoplot}}} +}} + diff --git a/man/serocalculator.Rd b/man/serocalculator.Rd index 33599f03..e30fa759 100644 --- a/man/serocalculator.Rd +++ b/man/serocalculator.Rd @@ -5,27 +5,10 @@ \alias{serocalculator-package} \title{Estimating Infection Rates from Serological Data} \description{ -This package translates antibody levels measured in a (cross-sectional) population sample into an -estimate of the frequency with which seroconversions (infections) occur in the sampled population. - -The API for this package includes the following functions:\tabular{ll}{ - Function Name \tab Purpose \cr - \code{\link[=load_pop_data]{load_pop_data()}} \tab loading cross-sectional antibody survey data \cr - \code{\link[=check_pop_data]{check_pop_data()}} \tab checking antibody data \cr - \code{\link[=summary.pop_data]{summary.pop_data()}} \tab numerical summaries of antibody data \cr - \code{\link[=autoplot.pop_data]{autoplot.pop_data()}} \tab graphs of antibody data distributions \cr - \code{\link[=load_curve_params]{load_curve_params()}} \tab loading antibody decay curve models \cr - \code{\link[=autoplot.curve_params]{autoplot.curve_params()}} \tab graphing antibody decay curves \cr - \code{\link[=llik]{llik()}} \tab computing log-likelihoods \cr - \code{\link[=graph.loglik]{graph.loglik()}} \tab graphing log-likelihood functions \cr - \code{\link[=autoplot.seroincidence]{autoplot.seroincidence()}} \tab graphing log-likelihood functions \cr - \code{\link[=autoplot.seroincidence.by]{autoplot.seroincidence.by()}} \tab graphing log-likelihood functions \cr - \code{\link[=est.incidence]{est.incidence()}} \tab estimating incidence rates \cr - \code{\link[=est.incidence.by]{est.incidence.by()}} \tab estimating incidence rates by strata \cr - \code{\link[=summary.seroincidence.by]{summary.seroincidence.by()}} \tab summarizing stratified incidence rate estimates \cr - \code{\link[=autoplot.summary.seroincidence.by]{autoplot.summary.seroincidence.by()}} \tab graphing incidence rate estimates \cr - \code{\link[=sim.cs]{sim.cs()}} \tab simulating cross-sectional population antibody data using longitudinal seroresponse models \cr -} +This package translates antibody levels measured in a +(cross-sectional) population sample into an +estimate of the frequency with which seroconversions (infections) +occur in the sampled population. } \details{ _PACKAGE @@ -33,48 +16,88 @@ _PACKAGE \references{ \emph{\strong{Methods for estimating seroincidence}} \itemize{ -\item Teunis, P. F. M., and J. C. H. van Eijkeren. "Estimation of seroconversion rates for infectious diseases: Effects of age and noise." Statistics in Medicine 39, no. 21 (2020): 2799-2814. -\item Teunis, P. F. M., J. C. H. van Eijkeren, W. F. de Graaf, A. Bonačić Marinović, and M. E. E. Kretzschmar. "Linking the seroresponse to infection to within-host heterogeneity in antibody production." Epidemics 16 (2016): 33-39. +\item Teunis, P. F. M., and J. C. H. van Eijkeren. +"Estimation of seroconversion rates for infectious diseases: +Effects of age and noise." +Statistics in Medicine 39, no. 21 (2020): 2799-2814. +\item Teunis, P. F. M., J. C. H. van Eijkeren, W. F. de Graaf, +A. Bonačić Marinović, and M. E. E. Kretzschmar. +"Linking the seroresponse to infection to within-host heterogeneity +in antibody production." Epidemics 16 (2016): 33-39. } \emph{\strong{Applications}} \itemize{ -\item Aiemjoy, Kristen, Jessica C. Seidman, Senjuti Saha, Sira Jam Munira, Mohammad Saiful Islam Sajib, Syed Muktadir Al Sium, Anik Sarkar et al. "Estimating typhoid incidence from community-based serosurveys: a multicohort study." The Lancet Microbe 3, no. 8 (2022): e578-e587. -\item Aiemjoy, Kristen, John Rumunu, Juma John Hassen, Kirsten E. Wiens, Denise Garrett, Polina Kamenskaya, Jason B. Harris et al. "Seroincidence of enteric fever, Juba, South Sudan." Emerging infectious diseases 28, no. 11 (2022): 2316. -\item Monge, S., Teunis, P. F., Friesema, I., Franz, E., Ang, W., van Pelt, W., Mughini-Gras, L. -"Immune response-eliciting exposure to Campylobacter vastly exceeds the incidence of clinically -overt campylobacteriosis but is associated with similar risk factors: A nationwide serosurvey in the Netherlands" +\item Aiemjoy, Kristen, Jessica C. Seidman, Senjuti Saha, +Sira Jam Munira, Mohammad Saiful Islam Sajib, +Syed Muktadir Al Sium, Anik Sarkar et al. +"Estimating typhoid incidence from community-based serosurveys: +a multicohort study." The Lancet Microbe 3, no. 8 (2022): e578-e587. +\item Aiemjoy, Kristen, John Rumunu, Juma John Hassen, +Kirsten E. Wiens, Denise Garrett, Polina Kamenskaya, +Jason B. Harris et al. "Seroincidence of enteric fever, Juba, +South Sudan." Emerging infectious diseases 28, no. 11 (2022): 2316. +\item Monge, S., Teunis, P. F., Friesema, I., Franz, E., Ang, W., +van Pelt, W., Mughini-Gras, L. +"Immune response-eliciting exposure to Campylobacter vastly exceeds +the incidence of clinically +overt campylobacteriosis but is associated with similar +risk factors: A nationwide serosurvey in the Netherlands" Journal of Infection, 2018, 1--7, doi:10.1016/j.jinf.2018.04.016 \item Kretzschmar, M., Teunis, P. F., Pebody, R. G. -"Incidence and reproduction numbers of pertussis: estimates from serological and social contact data in five European countries" -PLoS Medicine 7, no. 6 (June 1, 2010):e1000291. doi:10.1371/journal.pmed.1000291. -\item Simonsen, J., Strid, M. A., Molbak, K., Krogfelt, K. A., Linneberg, A., Teunis, P. -"Sero-epidemiology as a tool to study the incidence of Salmonella infections in humans" -Epidemiology and Infection 136, no. 7 (July 1, 2008): 895--902. doi:10.1017/S0950268807009314 -\item Simonsen, J., Teunis, P. F., van Pelt, W., van Duynhoven, Y., Krogfelt, K. A., Sadkowska-Todys, M., Molbak K. -"Usefulness of seroconversion rates for comparing infection pressures between countries" -Epidemiology and Infection, April 12, 2010, 1--8. doi:10.1017/S0950268810000750. -\item Falkenhorst, G., Simonsen, J., Ceper, T. H., van Pelt, W., de Valk, H., Sadkowska-Todys, M., Zota, L., Kuusi, M., Jernberg, C., Rota, M. C., van Duynhoven, Y. T., Teunis, P. F., Krogfelt, K. A., Molbak, K. -"Serological cross-sectional studies on salmonella incidence in eight European countries: no correlation with incidence of reported cases" -BMC Public Health 12, no. 1 (July 15, 2012): 523--23. doi:10.1186/1471-2458-12-523. -\item Teunis, P. F., Falkenhorst, G., Ang, C. W., Strid, M. A., De Valk, H., Sadkowska-Todys, M., Zota, L., Kuusi, M., Rota, M. C., Simonsen, J. B., Molbak, K., Van Duynhoven, Y. T., van Pelt, W. +"Incidence and reproduction numbers of pertussis: estimates from +serological and social contact data in five European countries" +PLoS Medicine 7, no. 6 (June 1, 2010):e1000291. +doi:10.1371/journal.pmed.1000291. +\item Simonsen, J., Strid, M. A., Molbak, K., Krogfelt, K. A., +Linneberg, A., Teunis, P. +"Sero-epidemiology as a tool to study the incidence of +Salmonella infections in humans" +Epidemiology and Infection 136, no. 7 (July 1, 2008): 895--902. +doi:10.1017/S0950268807009314 +\item Simonsen, J., Teunis, P. F., van Pelt, W., van Duynhoven, Y., +Krogfelt, K. A., Sadkowska-Todys, M., Molbak K. +"Usefulness of seroconversion rates for comparing infection +pressures between countries" +Epidemiology and Infection, April 12, 2010, 1--8. +doi:10.1017/S0950268810000750. +\item Falkenhorst, G., Simonsen, J., Ceper, T. H., van Pelt, W., +de Valk, H., Sadkowska-Todys, M., Zota, L., Kuusi, M., +Jernberg, C., Rota, M. C., van Duynhoven, Y. T., Teunis, P. F., +Krogfelt, K. A., Molbak, K. +"Serological cross-sectional studies on salmonella incidence in +eight European countries: no correlation with incidence of +reported cases" +BMC Public Health 12, no. 1 (July 15, 2012): 523--23. +doi:10.1186/1471-2458-12-523. +\item Teunis, P. F., Falkenhorst, G., Ang, C. W., Strid, M. A., +De Valk, H., Sadkowska-Todys, M., Zota, L., Kuusi, M., Rota, M. C., Simonsen, J. B., Molbak, K., Van Duynhoven, Y. T., van Pelt, W. "Campylobacter seroconversion rates in selected countries in the European Union" -Epidemiology and Infection 141, no. 10 (2013): 2051--57. doi:10.1017/S0950268812002774. -\item de Melker, H. E., Versteegh, F. G., Schellekens, J. F., Teunis, P. F., Kretzschmar, M. -"The incidence of Bordetella pertussis infections estimated in the population from a combination of serological surveys" +Epidemiology and Infection 141, no. 10 (2013): 2051--57. +doi:10.1017/S0950268812002774. +\item de Melker, H. E., Versteegh, F. G., Schellekens, J. F., +Teunis, P. F., Kretzschmar, M. +"The incidence of Bordetella pertussis infections estimated in the +population from a combination of serological surveys" The Journal of Infection 53, no. 2 (August 1, 2006): 106--13. doi:10.1016/j.jinf.2005.10.020 } \emph{\strong{Quantification of seroresponse}} \itemize{ \item de Graaf, W. F., Kretzschmar, M. E., Teunis, P. F., Diekmann, O. -"A two-phase within-host model for immune response and its application to serological profiles of pertussis" +"A two-phase within-host model for immune response and its +application to serological profiles of pertussis" Epidemics 9 (2014):1--7. doi:10.1016/j.epidem.2014.08.002. -\item Berbers, G. A., van de Wetering, M. S., van Gageldonk, P. G., Schellekens, J. F., Versteegh, F. G., Teunis, P.F. -"A novel method for evaluating natural and vaccine induced serological responses to Bordetella pertussis antigens" -Vaccine 31, no. 36 (August 12, 2013): 3732--38. doi:10.1016/j.vaccine.2013.05.073. -\item Versteegh, F. G., Mertens, P. L., de Melker, H. E., Roord, J. J., Schellekens, J. F., Teunis, P. F. -"Age-specific long-term course of IgG antibodies to pertussis toxin after symptomatic infection with Bordetella pertussis" +\item Berbers, G. A., van de Wetering, M. S., van Gageldonk, P. G., +Schellekens, J. F., Versteegh, F. G., Teunis, P.F. +"A novel method for evaluating natural and vaccine induced +serological responses to Bordetella pertussis antigens" +Vaccine 31, no. 36 (August 12, 2013): 3732--38. +doi:10.1016/j.vaccine.2013.05.073. +\item Versteegh, F. G., Mertens, P. L., de Melker, H. E., Roord, J. J., +Schellekens, J. F., Teunis, P. F. +"Age-specific long-term course of IgG antibodies to pertussis toxin +after symptomatic infection with Bordetella pertussis" Epidemiology and Infection 133, no. 4 (August 1, 2005): 737--48. \item Teunis, P. F., van der Heijden, O. G., de Melker, H. E., Schellekens, J. F., Versteegh, F. G., Kretzschmar, M. E. "Kinetics of the IgG antibody response to pertussis toxin after infection with B. pertussis" diff --git a/man/sim.cs.Rd b/man/sim.cs.Rd index 74b2ccd0..01bae0a6 100644 --- a/man/sim.cs.Rd +++ b/man/sim.cs.Rd @@ -45,14 +45,18 @@ sim.cs( \item{add.noise}{a \code{\link[=logical]{logical()}} indicating whether to add biological and measurement noise} -\item{curve_params}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters from the Bayesian posterior distribution of a longitudinal decay curve model. The parameter columns must be named: +\item{curve_params}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters +from the Bayesian posterior distribution of a longitudinal decay curve model. +The parameter columns must be named: \itemize{ -\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype combinations +\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype +combinations \item \code{iter}: an \code{\link[=integer]{integer()}} vector indicating MCMC sampling iterations \item \code{y0}: baseline antibody level at $t=0$ ($y(t=0)$) \item \code{y1}: antibody peak level (ELISA units) \item \code{t1}: duration of infection -\item \code{alpha}: antibody decay rate (1/days for the current longitudinal parameter sets) +\item \code{alpha}: antibody decay rate +(1/days for the current longitudinal parameter sets) \item \code{r}: shape factor of antibody decay }} diff --git a/man/sim.cs.multi.Rd b/man/sim.cs.multi.Rd index 3f8cd593..fe2204cc 100644 --- a/man/sim.cs.multi.Rd +++ b/man/sim.cs.multi.Rd @@ -53,14 +53,18 @@ sim.cs.multi( \item \code{"long"} (one measurement per row) or \item \code{"wide"} (one serum sample per row) }} - \item{\code{curve_params}}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters from the Bayesian posterior distribution of a longitudinal decay curve model. The parameter columns must be named: + \item{\code{curve_params}}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters +from the Bayesian posterior distribution of a longitudinal decay curve model. +The parameter columns must be named: \itemize{ -\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype combinations +\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype +combinations \item \code{iter}: an \code{\link[=integer]{integer()}} vector indicating MCMC sampling iterations \item \code{y0}: baseline antibody level at $t=0$ ($y(t=0)$) \item \code{y1}: antibody peak level (ELISA units) \item \code{t1}: duration of infection -\item \code{alpha}: antibody decay rate (1/days for the current longitudinal parameter sets) +\item \code{alpha}: antibody decay rate +(1/days for the current longitudinal parameter sets) \item \code{r}: shape factor of antibody decay }} }} diff --git a/man/stratify_data.Rd b/man/stratify_data.Rd index 7dacbcbb..fba7c6f0 100644 --- a/man/stratify_data.Rd +++ b/man/stratify_data.Rd @@ -15,20 +15,27 @@ stratify_data( ) } \arguments{ -\item{curve_params}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters from the Bayesian posterior distribution of a longitudinal decay curve model. The parameter columns must be named: +\item{curve_params}{a \code{\link[=data.frame]{data.frame()}} containing MCMC samples of parameters +from the Bayesian posterior distribution of a longitudinal decay curve model. +The parameter columns must be named: \itemize{ -\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype combinations +\item \code{antigen_iso}: a \code{\link[=character]{character()}} vector indicating antigen-isotype +combinations \item \code{iter}: an \code{\link[=integer]{integer()}} vector indicating MCMC sampling iterations \item \code{y0}: baseline antibody level at $t=0$ ($y(t=0)$) \item \code{y1}: antibody peak level (ELISA units) \item \code{t1}: duration of infection -\item \code{alpha}: antibody decay rate (1/days for the current longitudinal parameter sets) +\item \code{alpha}: antibody decay rate +(1/days for the current longitudinal parameter sets) \item \code{r}: shape factor of antibody decay }} -\item{noise_params}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) containing the following variables, specifying noise parameters for each antigen isotype: +\item{noise_params}{a \code{\link[=data.frame]{data.frame()}} (or \code{\link[tibble:tibble]{tibble::tibble()}}) +containing the following variables, +specifying noise parameters for each antigen isotype: \itemize{ -\item \code{antigen_iso}: antigen isotype whose noise parameters are being specified on each row +\item \code{antigen_iso}: antigen isotype whose noise parameters are being specified +on each row \item \code{nu}: biological noise \item \code{eps}: measurement noise \item \code{y.low}: lower limit of detection for the current antigen isotype diff --git a/tests/spelling.R b/tests/spelling.R index e55675da..6713838f 100644 --- a/tests/spelling.R +++ b/tests/spelling.R @@ -1,3 +1,3 @@ if(requireNamespace('spelling', quietly = TRUE)) - spelling::spell_check_test(vignettes = TRUE, error = TRUE, + spelling::spell_check_test(vignettes = TRUE, error = FALSE, skip_on_cran = TRUE) diff --git a/tests/testthat/_snaps/as_curve_params.md b/tests/testthat/_snaps/as_curve_params.md index 5d1b1bd5..f4057784 100644 --- a/tests/testthat/_snaps/as_curve_params.md +++ b/tests/testthat/_snaps/as_curve_params.md @@ -120,5 +120,6 @@ AnQxAAQACQAAAAVhbHBoYQAEAAkAAAABcgAABAIAAAABAAQACQAAAAxhbnRpZ2VuX2lzb3MA AAMNAAAABQAAAAEAAAACAAAAAwAAAAQAAAAFAAAEAgAAAf8AAAAQAAAABQAEAAkAAAAISGx5 RV9JZ0EABAAJAAAACEhseUVfSWdHAAQACQAAAAdMUFNfSWdBAAQACQAAAAdMUFNfSWdHAAQA - CQAAAAZWaV9JZ0cAAAQCAAAC/wAAABAAAAABAAQACQAAAAZmYWN0b3IAAAD+AAAA/g== + CQAAAAZWaV9JZ0cAAAQCAAAC/wAAABAAAAABAAQACQAAAAZmYWN0b3IAAAD+AAAEAgAAAAEA + BAAJAAAADWJpb21hcmtlcl92YXIAAAAQAAAAAQAEAAkAAAALYW50aWdlbl9pc28AAAD+ diff --git a/tests/testthat/_snaps/as_curve_params/curve-data.csv b/tests/testthat/_snaps/as_curve_params/curve-data.csv new file mode 100644 index 00000000..263775ae --- /dev/null +++ b/tests/testthat/_snaps/as_curve_params/curve-data.csv @@ -0,0 +1,101 @@ +antigen_iso,iter,y0,y1,t1,alpha,r +HlyE_IgA,1,2.48226,63.4814,9.51957,5.8075e-4,1.74504 +HlyE_IgG,1,3.03648,163.893,6.54931,0.00457307,1.17417 +LPS_IgA,1,0.748207,103.122,4.97771,0.00307839,1.58393 +LPS_IgG,1,0.940505,319.669,6.13807,0.00166412,1.41179 +Vi_IgG,1,8.4628,4347.99,3.07151,3.40357e-5,1.05744 +HlyE_IgA,2,3.85638,288.162,1.26542,4.59482e-4,2.66301 +HlyE_IgG,2,1.81736,153.629,10.7868,9.20577e-4,1.29748 +LPS_IgA,2,1.75926,851.907,2.49327,1.26241e-4,2.90784 +LPS_IgG,2,0.433893,20.5905,4.29147,0.00121652,1.36993 +Vi_IgG,2,18.8354,345.078,3.48079,1.42191e-4,1.02415 +HlyE_IgA,3,2.00641,431.819,5.46487,2.77272e-4,1.61055 +HlyE_IgG,3,4.93741,143.455,8.58692,0.00293646,1.13118 +LPS_IgA,3,1.5608,154.713,2.99669,6.18879e-4,2.38064 +LPS_IgG,3,5.68858,109.978,1.09117,1.02102e-4,2.67518 +Vi_IgG,3,27.739,1599.76,2.79969,1.5827e-4,1.03128 +HlyE_IgA,4,2.29369,30.5523,6.85193,0.00126803,1.87014 +HlyE_IgG,4,1.3733900000000001,53.3082,11.0627,1.85056e-4,2.12891 +LPS_IgA,4,1.43717,337.982,2.03,8.44129e-5,2.17232 +LPS_IgG,4,0.618952,100.897,5.61409,0.00432757,1.49847 +Vi_IgG,4,1.43201,510.143,13.1698,1.43093e-5,1.84891 +HlyE_IgA,5,0.816123,160.249,6.89268,0.00139768,1.39745 +HlyE_IgG,5,1.67354,557.058,7.36457,5.05441e-4,1.97494 +LPS_IgA,5,0.686836,31.0828,3.21353,3.76136e-4,1.68504 +LPS_IgG,5,0.647737,407.455,2.92829,7.60085e-4,1.33092 +Vi_IgG,5,0.764798,1416.55,22.2927,1.20457e-5,1.52949 +HlyE_IgA,6,4.50001,524.75,1.47307,2.94074e-4,2.26046 +HlyE_IgG,6,2.08975,957.184,23.9275,6.5467e-4,1.73952 +LPS_IgA,6,0.654862,39.1707,4.35167,8.0372e-4,1.44848 +LPS_IgG,6,1.19912,163.357,1.95468,5.01458e-4,1.22989 +Vi_IgG,6,6.82922,668.296,3.24755,4.40317e-5,1.07209 +HlyE_IgA,7,0.954787,30.7757,17.036,0.00459319,1.54839 +HlyE_IgG,7,3.65695,696.886,4.55038,8.87688e-4,1.28947 +LPS_IgA,7,3.56253,527.923,1.61473,3.04477e-5,3.19588 +LPS_IgG,7,8.52327,69.7724,0.68349,1.23209e-4,2.82054 +Vi_IgG,7,8.66238,490.685,1.85229,3.01779e-5,1.06247 +HlyE_IgA,8,2.18397,41.3351,3.97109,0.00233809,1.47712 +HlyE_IgG,8,1.56706,49.6963,7.28684,2.0987e-4,1.6087 +LPS_IgA,8,1.90888,579.924,2.53354,1.42503e-4,2.37684 +LPS_IgG,8,2.20698,92.9535,0.942893,7.36317e-4,2.46131 +Vi_IgG,8,48.2299,293.986,0.564755,3.62345e-5,1.0182 +HlyE_IgA,9,96.4399,248.453,0.497351,4.67003e-5,8.49638 +HlyE_IgG,9,4.35328,284.239,3.9641,0.00104004,1.22787 +LPS_IgA,9,1.90063,562.177,2.67798,2.35436e-4,2.24837 +LPS_IgG,9,0.533052,32.0628,3.66596,0.0037518,1.16723 +Vi_IgG,9,32.6384,1491.45,1.58666,8.46132e-5,1.03933 +HlyE_IgA,10,3.30934,319.258,1.28441,4.48299e-4,1.7493 +HlyE_IgG,10,2.06421,55.1973,7.96327,2.07444e-4,2.09051 +LPS_IgA,10,0.277746,4.06174,3.47415,0.0405528,1.52969 +LPS_IgG,10,1.51317,113.998,1.20596,6.72707e-4,1.77359 +Vi_IgG,10,277.372,388.961,0.259715,1.52336e-4,1.00258 +HlyE_IgA,11,0.372154,9.85591,11.1798,0.0018439,1.25153 +HlyE_IgG,11,1.21291,142.482,10.1788,5.41775e-4,1.64794 +LPS_IgA,11,0.625559,767.291,2.9934,0.00402444,2.05997 +LPS_IgG,11,9.98187,389.894,0.780147,6.88363e-5,2.41826 +Vi_IgG,11,14.8223,1563.58,2.42263,2.43812e-5,1.07932 +HlyE_IgA,12,3.04926,58.5691,1.2496,1.65952e-4,1.8784 +HlyE_IgG,12,2.65381,125.103,9.43966,0.00322191,1.1688 +LPS_IgA,12,2.05564,894.809,2.43653,6.32314e-5,2.71435 +LPS_IgG,12,15.2404,328.443,0.498514,3.90956e-5,2.67242 +Vi_IgG,12,7.2945,915.872,4.43754,3.42337e-5,1.13689 +HlyE_IgA,13,0.814143,23.5781,9.11408,8.01609e-4,1.27396 +HlyE_IgG,13,1.7148,154.135,13.0894,0.00117192,1.42242 +LPS_IgA,13,0.834687,8.76744,3.06246,0.00290858,1.86291 +LPS_IgG,13,0.581101,54.0885,4.68085,0.00220011,1.29866 +Vi_IgG,13,1.88268,503.456,7.94768,3.74989e-5,1.71411 +HlyE_IgA,14,1.7329,195.272,4.46496,0.00212248,1.60775 +HlyE_IgG,14,3.02769,512.686,8.49536,9.76984e-4,1.22367 +LPS_IgA,14,3.45558,113.417,2.16853,9.64496e-5,2.69646 +LPS_IgG,14,1.31537,140.974,2.17098,0.00134016,1.46963 +Vi_IgG,14,10.4619,520.314,1.73096,6.00634e-5,1.0501 +HlyE_IgA,15,0.324205,15.2392,8.14611,0.00206185,1.13545 +HlyE_IgG,15,2.19975,367.676,16.037,3.8458e-4,1.45282 +LPS_IgA,15,1.24141,6.95947,2.56553,0.00138079,1.89075 +LPS_IgG,15,14.3329,365.608,0.446341,4.34935e-5,4.82298 +Vi_IgG,15,2.97466,529.611,13.2801,1.36865e-5,1.24711 +HlyE_IgA,16,3.91941,138.598,5.73605,5.20896e-4,1.83739 +HlyE_IgG,16,1.29501,94.0545,5.45728,0.00100652,1.38485 +LPS_IgA,16,1.68487,114.235,2.49778,4.10202e-4,2.33134 +LPS_IgG,16,0.949917,107.109,2.11729,5.59124e-4,1.52846 +Vi_IgG,16,2.19421,727.055,8.82501,3.07599e-5,1.39757 +HlyE_IgA,17,2.96667,17.7003,12.1549,0.00110521,1.33446 +HlyE_IgG,17,1.37419,420.52,27.3512,1.83195e-4,2.2181 +LPS_IgA,17,2.95752,4736.5,2.14502,1.1537e-4,2.69753 +LPS_IgG,17,3.78808,378.497,1.03972,2.59878e-4,1.51897 +Vi_IgG,17,23.3291,122.167,2.55241,4.9078e-5,1.0338 +HlyE_IgA,18,12.6872,34.6298,4.87057,5.03901e-4,2.87854 +HlyE_IgG,18,1.14206,117.117,9.70208,2.28263e-4,2.95589 +LPS_IgA,18,2.59052,499.625,2.04702,3.9643e-4,2.22583 +LPS_IgG,18,1.14575,219.817,2.23021,4.63748e-4,1.58808 +Vi_IgG,18,8.78688,677.992,5.29091,1.63194e-4,1.22066 +HlyE_IgA,19,0.633407,670.782,10.7817,0.00194229,1.37775 +HlyE_IgG,19,0.557342,139.745,21.1086,1.71501e-4,1.7877 +LPS_IgA,19,1.17559,2652.75,3.30762,7.89449e-5,2.28541 +LPS_IgG,19,1.27842,68.5481,2.53069,6.06813e-4,1.49596 +Vi_IgG,19,25.9611,552.645,1.46112,9.88853e-5,1.02806 +HlyE_IgA,20,15.9786,350.389,3.57177,3.69696e-4,3.42523 +HlyE_IgG,20,3.85276,280.773,4.45087,0.00104854,1.21638 +LPS_IgA,20,2.00689,8.14268,2.43644,5.09893e-4,2.04843 +LPS_IgG,20,0.746152,123.385,3.067,9.18416e-4,1.39196 +Vi_IgG,20,0.827451,714.012,11.481,1.37235e-5,1.97522 diff --git a/tests/testthat/_snaps/class_attributes.md b/tests/testthat/_snaps/class_attributes.md new file mode 100644 index 00000000..47b50a72 --- /dev/null +++ b/tests/testthat/_snaps/class_attributes.md @@ -0,0 +1,430 @@ +# `get_id()` works + + c("B1", "B10", "B100", "B101", "B102", "B103", "B104", "B105", + "B106", "B107", "B108", "B109", "B11", "B110", "B111", "B112", + "B113", "B114", "B115", "B116", "B117", "B118", "B119", "B12", + "B120", "B121", "B122", "B123", "B124", "B125", "B126", "B127", + "B128", "B129", "B13", "B130", "B131", "B132", "B133", "B134", + "B135", "B136", "B137", "B138", "B139", "B14", "B140", "B141", + "B142", "B143", "B144", "B145", "B146", "B147", "B148", "B149", + "B15", "B150", "B151", "B152", "B153", "B154", "B155", "B156", + "B157", "B158", "B159", "B16", "B160", "B161", "B162", "B163", + "B164", "B165", "B166", "B167", "B168", "B169", "B17", "B170", + "B171", "B172", "B173", "B174", "B175", "B176", "B177", "B178", + "B179", "B18", "B180", "B181", "B182", "B183", "B184", "B185", + "B186", "B187", "B188", "B189", "B19", "B190", "B191", "B192", + "B193", "B194", "B195", "B196", "B197", "B198", "B199", "B2", + "B20", "B200", "B201", "B202", "B203", "B204", "B205", "B206", + "B207", "B208", "B209", "B21", "B210", "B211", "B212", "B213", + "B214", "B215", "B216", "B217", "B218", "B219", "B22", "B220", + "B221", "B222", "B223", "B224", "B225", "B226", "B227", "B228", + "B229", "B23", "B230", "B231", "B232", "B233", "B234", "B235", + "B236", "B237", "B238", "B239", "B24", "B240", "B241", "B242", + "B243", "B244", "B245", "B246", "B247", "B248", "B249", "B25", + "B250", "B251", "B252", "B253", "B254", "B255", "B256", "B257", + "B258", "B259", "B26", "B260", "B261", "B262", "B263", "B264", + "B265", "B266", "B267", "B268", "B269", "B27", "B270", "B271", + "B272", "B273", "B274", "B275", "B276", "B277", "B278", "B279", + "B28", "B280", "B281", "B282", "B283", "B284", "B285", "B286", + "B287", "B288", "B289", "B29", "B290", "B291", "B292", "B293", + "B294", "B295", "B296", "B297", "B298", "B299", "B3", "B30", + "B300", "B301", "B302", "B303", "B304", "B305", "B306", "B307", + "B308", "B309", "B31", "B310", "B311", "B312", "B313", "B314", + "B315", "B316", "B317", "B318", "B319", "B32", "B320", "B321", + "B322", "B323", "B324", "B325", "B326", "B327", "B328", "B329", + "B33", "B330", "B331", "B332", "B333", "B334", "B335", "B336", + "B337", "B338", "B339", "B34", "B340", "B341", "B342", "B343", + "B344", "B345", "B346", "B347", "B348", "B349", "B35", "B350", + "B351", "B352", "B353", "B354", "B355", "B356", "B357", "B358", + "B359", "B36", "B360", "B361", "B362", "B363", "B364", "B365", + "B366", "B367", "B368", "B369", "B37", "B370", "B371", "B372", + "B373", "B374", "B375", "B376", "B377", "B378", "B379", "B38", + "B380", "B381", "B382", "B383", "B384", "B385", "B386", "B387", + "B388", "B389", "B39", "B390", "B391", "B392", "B393", "B394", + "B395", "B396", "B397", "B398", "B399", "B4", "B40", "B400", + "B401", "B402", "B403", "B404", "B405", "B406", "B407", "B408", + "B409", "B41", "B410", "B411", "B412", "B413", "B414", "B415", + "B416", "B417", "B418", "B419", "B42", "B420", "B421", "B422", + "B423", "B424", "B425", "B426", "B427", "B428", "B429", "B43", + "B430", "B431", "B432", "B433", "B434", "B435", "B436", "B437", + "B438", "B439", "B44", "B440", "B441", "B442", "B443", "B444", + "B445", "B446", "B447", "B448", "B449", "B45", "B450", "B451", + "B452", "B453", "B454", "B455", "B456", "B457", "B458", "B459", + "B46", "B460", "B461", "B462", "B463", "B464", "B465", "B466", + "B467", "B468", "B469", "B47", "B470", "B471", "B472", "B473", + "B474", "B475", "B476", "B477", "B478", "B479", "B48", "B480", + "B481", "B482", "B483", "B484", "B485", "B486", "B487", "B488", + "B489", "B49", "B490", "B491", "B492", "B493", "B494", "B495", + "B496", "B497", "B498", "B499", "B5", "B50", "B500", "B501", + "B502", "B503", "B504", "B505", "B506", "B507", "B508", "B509", + "B51", "B510", "B511", "B512", "B513", "B514", "B515", "B516", + "B517", "B518", "B519", "B52", "B520", "B521", "B522", "B523", + "B524", "B525", "B526", "B527", "B528", "B529", "B53", "B530", + "B531", "B532", "B533", "B534", "B535", "B536", "B537", "B538", + "B539", "B54", "B540", "B541", "B542", "B543", "B544", "B545", + "B546", "B547", "B548", "B549", "B55", "B550", "B551", "B552", + "B553", "B554", "B555", "B556", "B557", "B558", "B559", "B56", + "B560", "B561", "B562", "B563", "B564", "B565", "B566", "B567", + "B568", "B569", "B57", "B570", "B571", "B572", "B573", "B574", + "B575", "B576", "B577", "B578", "B579", "B58", "B580", "B581", + "B582", "B583", "B584", "B585", "B586", "B587", "B588", "B589", + "B59", "B590", "B591", "B592", "B593", "B594", "B595", "B596", + "B597", "B598", "B599", "B6", "B60", "B600", "B601", "B602", + "B603", "B604", "B605", "B606", "B607", "B608", "B609", "B61", + "B610", "B611", "B612", "B613", "B614", "B615", "B616", "B617", + "B618", "B619", "B62", "B620", "B621", "B622", "B623", "B624", + "B625", "B626", "B627", "B628", "B629", "B63", "B630", "B631", + "B632", "B633", "B634", "B635", "B636", "B637", "B638", "B639", + "B64", "B640", "B641", "B642", "B643", "B644", "B645", "B646", + "B647", "B648", "B649", "B65", "B650", "B651", "B652", "B653", + "B654", "B655", "B656", "B657", "B658", "B659", "B66", "B660", + "B661", "B662", "B663", "B664", "B665", "B666", "B667", "B668", + "B669", "B67", "B670", "B671", "B672", "B673", "B674", "B675", + "B676", "B677", "B678", "B679", "B68", "B680", "B681", "B682", + "B683", "B684", "B685", "B686", "B687", "B688", "B689", "B69", + "B690", "B691", "B692", "B693", "B694", "B695", "B696", "B697", + "B698", "B699", "B7", "B70", "B700", "B701", "B702", "B703", + "B704", "B705", "B706", "B707", "B708", "B709", "B71", "B710", + "B711", "B712", "B713", "B714", "B715", "B716", "B717", "B718", + "B719", "B72", "B720", "B721", "B722", "B723", "B724", "B725", + "B726", "B727", "B728", "B729", "B73", "B730", "B731", "B732", + "B733", "B734", "B735", "B736", "B737", "B738", "B739", "B74", + "B740", "B741", "B742", "B743", "B744", "B745", "B746", "B747", + "B748", "B749", "B75", "B750", "B751", "B752", "B753", "B754", + "B755", "B756", "B757", "B758", "B759", "B76", "B760", "B761", + "B762", "B763", "B764", "B765", "B766", "B767", "B768", "B769", + "B77", "B770", "B771", "B772", "B773", "B774", "B775", "B776", + "B777", "B778", "B779", "B78", "B780", "B781", "B782", "B783", + "B784", "B785", "B786", "B787", "B788", "B789", "B79", "B790", + "B791", "B792", "B793", "B794", "B795", "B796", "B797", "B798", + "B799", "B8", "B80", "B800", "B801", "B802", "B81", "B82", "B83", + "B84", "B85", "B86", "B87", "B88", "B89", "B9", "B90", "B91", + "B92", "B93", "B94", "B95", "B96", "B97", "B98", "B99", "N1", + "N10", "N100", "N1000", "N1001", "N1002", "N1003", "N1004", "N1005", + "N1006", "N1007", "N1008", "N1009", "N101", "N1010", "N1011", + "N1012", "N1013", "N1014", "N1015", "N1016", "N1017", "N1018", + "N1019", "N102", "N1020", "N1021", "N1022", "N1023", "N1024", + "N1025", "N1026", "N1027", "N1028", "N1029", "N103", "N1030", + "N1031", "N1032", "N1033", "N1034", "N1035", "N1036", "N1037", + "N1038", "N1039", "N104", "N1040", "N1041", "N1042", "N1043", + "N1044", "N1045", "N1046", "N1047", "N1048", "N1049", "N105", + "N1050", "N1051", "N1052", "N1053", "N1054", "N1055", "N1056", + "N1057", "N1058", "N1059", "N106", "N1060", "N1061", "N1062", + "N1063", "N1064", "N1065", "N1066", "N1067", "N1068", "N1069", + "N107", "N1070", "N1071", "N1072", "N1073", "N1074", "N1075", + "N1076", "N1077", "N1078", "N1079", "N108", "N1080", "N1081", + "N1082", "N1083", "N1084", "N1085", "N1086", "N1087", "N1088", + "N1089", "N109", "N1090", "N1091", "N1092", "N1093", "N1094", + "N1095", "N1096", "N1097", "N1098", "N1099", "N11", "N110", "N1100", + "N1101", "N1102", "N1103", "N1104", "N1105", "N1106", "N1107", + "N1108", "N1109", "N111", "N1110", "N1111", "N1112", "N1113", + "N1114", "N1115", "N1116", "N1117", "N1118", "N1119", "N112", + "N1120", "N1121", "N1122", "N1123", "N1124", "N1125", "N1126", + "N1127", "N1128", "N1129", "N113", "N1130", "N1131", "N1132", + "N1133", "N1134", "N1135", "N1136", "N1137", "N1138", "N1139", + "N114", "N1140", "N1141", "N1142", "N1143", "N1144", "N1145", + "N1146", "N1147", "N1148", "N1149", "N115", "N1150", "N1151", + "N1152", "N1153", "N1154", "N1155", "N1156", "N1157", "N1158", + "N1159", "N116", "N1160", "N1161", "N1162", "N1163", "N1164", + "N1165", "N1166", "N1167", "N1168", "N1169", "N117", "N1170", + "N1171", "N1172", "N1173", "N1174", "N1175", "N1176", "N1177", + "N1178", "N1179", "N118", "N1180", "N1181", "N1182", "N1183", + "N1184", "N1185", "N1186", "N1187", "N1188", "N1189", "N119", + "N1190", "N1191", "N1192", "N1193", "N1194", "N1195", "N1196", + "N1197", "N1198", "N1199", "N12", "N120", "N1200", "N1201", "N1202", + "N1203", "N1204", "N1205", "N1206", "N1207", "N1208", "N1209", + "N121", "N1210", "N1211", "N1212", "N1213", "N1214", "N1215", + "N1216", "N1217", "N1218", "N1219", "N122", "N1220", "N1221", + "N1222", "N1223", "N1224", "N1225", "N1226", "N1227", "N1228", + "N1229", "N123", "N1230", "N1231", "N1232", "N1233", "N1234", + "N1235", "N1236", "N1237", "N1238", "N1239", "N124", "N1240", + "N1241", "N1242", "N1243", "N1244", "N1245", "N1246", "N1247", + "N1248", "N1249", "N125", "N1250", "N1251", "N1252", "N1253", + "N1254", "N1255", "N1256", "N1257", "N1258", "N1259", "N126", + "N1260", "N1261", "N1262", "N1263", "N1264", "N1265", "N1266", + "N1267", "N1268", "N1269", "N127", "N1270", "N1271", "N1272", + "N1273", "N1274", "N1275", "N1276", "N1277", "N1278", "N1279", + "N128", "N1280", "N1281", "N1282", "N1283", "N1284", "N1285", + "N1286", "N1287", "N1288", "N1289", "N129", "N1290", "N1291", + "N1292", "N1293", "N1294", "N1295", "N1296", "N1297", "N1298", + "N1299", "N13", "N130", "N1300", "N1301", "N1302", "N1303", "N1304", + "N1305", "N1306", "N1307", "N1308", "N1309", "N131", "N1310", + "N1311", "N1312", "N1313", "N1314", "N1315", "N1316", "N1317", + "N1318", "N1319", "N132", "N1320", "N1321", "N1322", "N1323", + "N1324", "N1325", "N1326", "N1327", "N1328", "N1329", "N133", + "N1330", "N1331", "N1332", "N1333", "N1334", "N1335", "N1336", + "N1337", "N1338", "N1339", "N134", "N1340", "N1341", "N1342", + "N1343", "N1344", "N1345", "N1346", "N1347", "N1348", "N1349", + "N135", "N1350", "N1351", "N1352", "N1353", "N1354", "N1355", + "N1356", "N1357", "N1358", "N1359", "N136", "N1360", "N1361", + "N1362", "N1363", "N1364", "N1365", "N1366", "N1367", "N1368", + "N1369", "N137", "N1370", "N1371", "N1372", "N1373", "N1374", + "N1375", "N1376", "N1377", "N1378", "N1379", "N138", "N1380", + "N1381", "N1382", "N1383", "N1384", "N1385", "N1386", "N1387", + "N1388", "N1389", "N139", "N1390", "N1391", "N1392", "N1393", + "N1394", "N1395", "N1396", "N1397", "N1398", "N1399", "N14", + "N140", "N1400", "N1401", "N1402", "N1403", "N1404", "N1405", + "N1406", "N1407", "N1408", "N1409", "N141", "N1410", "N1411", + "N1412", "N1413", "N1414", "N1415", "N1416", "N1417", "N1418", + "N1419", "N142", "N1420", "N1421", "N1422", "N1423", "N1424", + "N1425", "N1426", "N1427", "N1428", "N1429", "N143", "N1430", + "N1431", "N1432", "N1433", "N1434", "N1435", "N1436", "N1437", + "N1438", "N1439", "N144", "N1440", "N1441", "N1442", "N1443", + "N1444", "N1445", "N1446", "N1447", "N1448", "N1449", "N145", + "N1450", "N1451", "N1452", "N1453", "N1454", "N1455", "N1456", + "N1457", "N1458", "N1459", "N146", "N1460", "N1461", "N1462", + "N1463", "N1464", "N1465", "N1466", "N1467", "N1468", "N1469", + "N147", "N1470", "N1471", "N1472", "N1473", "N1474", "N1475", + "N1476", "N1477", "N1478", "N1479", "N148", "N1480", "N1481", + "N1482", "N1483", "N1484", "N1485", "N1486", "N1487", "N1488", + "N1489", "N149", "N1490", "N1491", "N1492", "N1493", "N1494", + "N1495", "N1496", "N1497", "N1498", "N1499", "N15", "N150", "N1500", + "N1501", "N1502", "N1503", "N1504", "N1505", "N1506", "N1507", + "N1508", "N1509", "N151", "N1510", "N1511", "N1512", "N1513", + "N1514", "N1515", "N1516", "N1517", "N1518", "N1519", "N152", + "N1520", "N1521", "N1522", "N1523", "N1524", "N1525", "N1526", + "N1527", "N1528", "N1529", "N153", "N1530", "N1531", "N1532", + "N1533", "N1534", "N1535", "N1536", "N1537", "N1538", "N1539", + "N154", "N1540", "N1541", "N1542", "N1543", "N1544", "N1545", + "N1546", "N155", "N156", "N157", "N158", "N159", "N16", "N160", + "N161", "N162", "N163", "N164", "N165", "N166", "N167", "N168", + "N169", "N17", "N170", "N171", "N172", "N173", "N174", "N175", + "N176", "N177", "N178", "N179", "N18", "N180", "N181", "N182", + "N183", "N184", "N185", "N186", "N187", "N188", "N189", "N19", + "N190", "N191", "N192", "N193", "N194", "N195", "N196", "N197", + "N198", "N199", "N2", "N20", "N200", "N201", "N202", "N203", + "N204", "N205", "N206", "N207", "N208", "N209", "N21", "N210", + "N211", "N212", "N213", "N214", "N215", "N216", "N217", "N218", + "N219", "N22", "N220", "N221", "N222", "N223", "N224", "N225", + "N226", "N227", "N228", "N229", "N23", "N230", "N231", "N232", + "N233", "N234", "N235", "N236", "N237", "N238", "N239", "N24", + "N240", "N241", "N242", "N243", "N244", "N245", "N246", "N247", + "N248", "N249", "N25", "N250", "N251", "N252", "N253", "N254", + "N255", "N256", "N257", "N258", "N259", "N26", "N260", "N261", + "N262", "N263", "N264", "N265", "N266", "N267", "N268", "N269", + "N27", "N270", "N271", "N272", "N273", "N274", "N275", "N276", + "N277", "N278", "N279", "N28", "N280", "N281", "N282", "N283", + "N284", "N285", "N286", "N287", "N288", "N289", "N29", "N290", + "N291", "N292", "N293", "N294", "N295", "N296", "N297", "N298", + "N299", "N3", "N30", "N300", "N301", "N302", "N303", "N304", + "N305", "N306", "N307", "N308", "N309", "N31", "N310", "N311", + "N312", "N313", "N314", "N315", "N316", "N317", "N318", "N319", + "N32", "N320", "N321", "N322", "N323", "N324", "N325", "N326", + "N327", "N328", "N329", "N33", "N330", "N331", "N332", "N333", + "N334", "N335", "N336", "N337", "N338", "N339", "N34", "N340", + "N341", "N342", "N343", "N344", "N345", "N346", "N347", "N348", + "N349", "N35", "N350", "N351", "N352", "N353", "N354", "N355", + "N356", "N357", "N358", "N359", "N36", "N360", "N361", "N362", + "N363", "N364", "N365", "N366", "N367", "N368", "N369", "N37", + "N370", "N371", "N372", "N373", "N374", "N375", "N376", "N377", + "N378", "N379", "N38", "N380", "N381", "N382", "N383", "N384", + "N385", "N386", "N387", "N388", "N389", "N39", "N390", "N391", + "N392", "N393", "N394", "N395", "N396", "N397", "N398", "N399", + "N4", "N40", "N400", "N401", "N402", "N403", "N404", "N405", + "N406", "N407", "N408", "N409", "N41", "N410", "N411", "N412", + "N413", "N414", "N415", "N416", "N417", "N418", "N419", "N42", + "N420", "N421", "N422", "N423", "N424", "N425", "N426", "N427", + "N428", "N429", "N43", "N430", "N431", "N432", "N433", "N434", + "N435", "N436", "N437", "N438", "N439", "N44", "N440", "N441", + "N442", "N443", "N444", "N445", "N446", "N447", "N448", "N449", + "N45", "N450", "N451", "N452", "N453", "N454", "N455", "N456", + "N457", "N458", "N459", "N46", "N460", "N461", "N462", "N463", + "N464", "N465", "N466", "N467", "N468", "N469", "N47", "N470", + "N471", "N472", "N473", "N474", "N475", "N476", "N477", "N478", + "N479", "N48", "N480", "N481", "N482", "N483", "N484", "N485", + "N486", "N487", "N488", "N489", "N49", "N490", "N491", "N492", + "N493", "N494", "N495", "N496", "N497", "N498", "N499", "N5", + "N50", "N500", "N501", "N502", "N503", "N504", "N505", "N506", + "N507", "N508", "N509", "N51", "N510", "N511", "N512", "N513", + "N514", "N515", "N516", "N517", "N518", "N519", "N52", "N520", + "N521", "N522", "N523", "N524", "N525", "N526", "N527", "N528", + "N529", "N53", "N530", "N531", "N532", "N533", "N534", "N535", + "N536", "N537", "N538", "N539", "N54", "N540", "N541", "N542", + "N543", "N544", "N545", "N546", "N547", "N548", "N549", "N55", + "N550", "N551", "N552", "N553", "N554", "N555", "N556", "N557", + "N558", "N559", "N56", "N560", "N561", "N562", "N563", "N564", + "N565", "N566", "N567", "N568", "N569", "N57", "N570", "N571", + "N572", "N573", "N574", "N575", "N576", "N577", "N578", "N579", + "N58", "N580", "N581", "N582", "N583", "N584", "N585", "N586", + "N587", "N588", "N589", "N59", "N590", "N591", "N592", "N593", + "N594", "N595", "N596", "N597", "N598", "N599", "N6", "N60", + "N600", "N601", "N602", "N603", "N604", "N605", "N606", "N607", + "N608", "N609", "N61", "N610", "N611", "N612", "N613", "N614", + "N615", "N616", "N617", "N618", "N619", "N62", "N620", "N621", + "N622", "N623", "N624", "N625", "N626", "N627", "N628", "N629", + "N63", "N630", "N631", "N632", "N633", "N634", "N635", "N636", + "N637", "N638", "N639", "N64", "N640", "N641", "N642", "N643", + "N644", "N645", "N646", "N647", "N648", "N649", "N65", "N650", + "N651", "N652", "N653", "N654", "N655", "N656", "N657", "N658", + "N659", "N66", "N660", "N661", "N662", "N663", "N664", "N665", + "N666", "N667", "N668", "N669", "N67", "N670", "N671", "N672", + "N673", "N674", "N675", "N676", "N677", "N678", "N679", "N68", + "N680", "N681", "N682", "N683", "N684", "N685", "N686", "N687", + "N688", "N689", "N69", "N690", "N691", "N692", "N693", "N694", + "N695", "N696", "N697", "N698", "N699", "N7", "N70", "N700", + "N701", "N702", "N703", "N704", "N705", "N706", "N707", "N708", + "N709", "N71", "N710", "N711", "N712", "N713", "N714", "N715", + "N716", "N717", "N718", "N719", "N72", "N720", "N721", "N722", + "N723", "N724", "N725", "N726", "N727", "N728", "N729", "N73", + "N730", "N731", "N732", "N733", "N734", "N735", "N736", "N737", + "N738", "N739", "N74", "N740", "N741", "N742", "N743", "N744", + "N745", "N746", "N747", "N748", "N749", "N75", "N750", "N751", + "N752", "N753", "N754", "N755", "N756", "N757", "N758", "N759", + "N76", "N760", "N761", "N762", "N763", "N764", "N765", "N766", + "N767", "N768", "N769", "N77", "N770", "N771", "N772", "N773", + "N774", "N775", "N776", "N777", "N778", "N779", "N78", "N780", + "N781", "N782", "N783", "N784", "N785", "N786", "N787", "N788", + "N789", "N79", "N790", "N791", "N792", "N793", "N794", "N795", + "N796", "N797", "N798", "N799", "N8", "N80", "N800", "N801", + "N802", "N803", "N804", "N805", "N806", "N807", "N808", "N809", + "N81", "N810", "N811", "N812", "N813", "N814", "N815", "N816", + "N817", "N818", "N819", "N82", "N820", "N821", "N822", "N823", + "N824", "N825", "N826", "N827", "N828", "N829", "N83", "N830", + "N831", "N832", "N833", "N834", "N835", "N836", "N837", "N838", + "N839", "N84", "N840", "N841", "N842", "N843", "N844", "N845", + "N846", "N847", "N848", "N849", "N85", "N850", "N851", "N852", + "N853", "N854", "N855", "N856", "N857", "N858", "N859", "N86", + "N860", "N861", "N862", "N863", "N864", "N865", "N866", "N867", + "N868", "N869", "N87", "N870", "N871", "N872", "N873", "N874", + "N875", "N876", "N877", "N878", "N879", "N88", "N880", "N881", + "N882", "N883", "N884", "N885", "N886", "N887", "N888", "N889", + "N89", "N890", "N891", "N892", "N893", "N894", "N895", "N896", + "N897", "N898", "N899", "N9", "N90", "N900", "N901", "N902", + "N903", "N904", "N905", "N906", "N907", "N908", "N909", "N91", + "N910", "N911", "N912", "N913", "N914", "N915", "N916", "N917", + "N918", "N919", "N92", "N920", "N921", "N922", "N923", "N924", + "N925", "N926", "N927", "N928", "N929", "N93", "N930", "N931", + "N932", "N933", "N934", "N935", "N936", "N937", "N938", "N939", + "N94", "N940", "N941", "N942", "N943", "N944", "N945", "N946", + "N947", "N948", "N949", "N95", "N950", "N951", "N952", "N953", + "N954", "N955", "N956", "N957", "N958", "N959", "N96", "N960", + "N961", "N962", "N963", "N964", "N965", "N966", "N967", "N968", + "N969", "N97", "N970", "N971", "N972", "N973", "N974", "N975", + "N976", "N977", "N978", "N979", "N98", "N980", "N981", "N982", + "N983", "N984", "N985", "N986", "N987", "N988", "N989", "N99", + "N990", "N991", "N992", "N993", "N994", "N995", "N996", "N997", + "N998", "N999", "P1", "P10", "P100", "P101", "P102", "P103", + "P104", "P105", "P106", "P107", "P108", "P109", "P11", "P110", + "P111", "P112", "P113", "P114", "P115", "P116", "P117", "P118", + "P119", "P12", "P120", "P121", "P122", "P123", "P124", "P125", + "P126", "P127", "P128", "P129", "P13", "P130", "P131", "P132", + "P133", "P134", "P135", "P136", "P137", "P138", "P139", "P14", + "P140", "P141", "P142", "P143", "P144", "P145", "P146", "P147", + "P148", "P149", "P15", "P150", "P151", "P152", "P153", "P154", + "P155", "P156", "P157", "P158", "P159", "P16", "P160", "P161", + "P162", "P163", "P164", "P165", "P166", "P167", "P168", "P169", + "P17", "P170", "P171", "P172", "P173", "P174", "P175", "P176", + "P177", "P178", "P179", "P18", "P180", "P181", "P182", "P183", + "P184", "P185", "P186", "P187", "P188", "P189", "P19", "P190", + "P191", "P192", "P193", "P194", "P195", "P196", "P197", "P198", + "P199", "P2", "P20", "P200", "P201", "P202", "P203", "P204", + "P205", "P206", "P207", "P208", "P209", "P21", "P210", "P211", + "P212", "P213", "P214", "P215", "P216", "P217", "P218", "P219", + "P22", "P220", "P221", "P222", "P223", "P224", "P225", "P226", + "P227", "P228", "P229", "P23", "P230", "P231", "P232", "P233", + "P234", "P235", "P236", "P237", "P238", "P239", "P24", "P240", + "P241", "P242", "P243", "P244", "P245", "P246", "P247", "P248", + "P249", "P25", "P250", "P251", "P252", "P253", "P254", "P255", + "P256", "P257", "P258", "P259", "P26", "P260", "P261", "P262", + "P263", "P264", "P265", "P266", "P267", "P268", "P269", "P27", + "P270", "P271", "P272", "P273", "P274", "P275", "P276", "P277", + "P278", "P279", "P28", "P280", "P281", "P282", "P283", "P284", + "P285", "P286", "P287", "P288", "P289", "P29", "P290", "P291", + "P292", "P293", "P294", "P295", "P296", "P297", "P298", "P299", + "P3", "P30", "P300", "P301", "P302", "P303", "P304", "P305", + "P306", "P307", "P308", "P309", "P31", "P310", "P311", "P312", + "P313", "P314", "P315", "P316", "P317", "P318", "P319", "P32", + "P320", "P321", "P322", "P323", "P324", "P325", "P326", "P327", + "P328", "P329", "P33", "P330", "P331", "P332", "P333", "P334", + "P335", "P336", "P337", "P338", "P339", "P34", "P340", "P341", + "P342", "P343", "P344", "P345", "P346", "P347", "P348", "P349", + "P35", "P350", "P351", "P352", "P353", "P354", "P355", "P356", + "P357", "P358", "P359", "P36", "P360", "P361", "P362", "P363", + "P364", "P365", "P366", "P367", "P368", "P369", "P37", "P370", + "P371", "P372", "P373", "P374", "P375", "P376", "P377", "P378", + "P379", "P38", "P380", "P381", "P382", "P383", "P384", "P385", + "P386", "P387", "P388", "P389", "P39", "P390", "P391", "P392", + "P393", "P394", "P395", "P396", "P397", "P398", "P399", "P4", + "P40", "P400", "P401", "P402", "P403", "P404", "P405", "P406", + "P407", "P408", "P409", "P41", "P410", "P411", "P412", "P413", + "P414", "P415", "P416", "P417", "P418", "P419", "P42", "P420", + "P421", "P422", "P423", "P424", "P425", "P426", "P427", "P428", + "P429", "P43", "P430", "P431", "P432", "P433", "P434", "P435", + "P436", "P437", "P438", "P439", "P44", "P440", "P441", "P442", + "P443", "P444", "P445", "P446", "P447", "P448", "P449", "P45", + "P450", "P451", "P452", "P453", "P454", "P455", "P456", "P457", + "P458", "P459", "P46", "P460", "P461", "P462", "P463", "P464", + "P465", "P466", "P467", "P468", "P469", "P47", "P470", "P471", + "P472", "P473", "P474", "P475", "P476", "P477", "P478", "P479", + "P48", "P480", "P481", "P482", "P483", "P484", "P485", "P486", + "P487", "P488", "P489", "P49", "P490", "P491", "P492", "P493", + "P494", "P495", "P496", "P497", "P498", "P499", "P5", "P50", + "P500", "P501", "P502", "P503", "P504", "P505", "P506", "P507", + "P508", "P509", "P51", "P510", "P511", "P512", "P513", "P514", + "P515", "P516", "P517", "P518", "P519", "P52", "P520", "P521", + "P522", "P523", "P524", "P525", "P526", "P527", "P528", "P529", + "P53", "P530", "P531", "P532", "P533", "P534", "P535", "P536", + "P537", "P538", "P539", "P54", "P540", "P541", "P542", "P543", + "P544", "P545", "P546", "P547", "P548", "P549", "P55", "P550", + "P551", "P552", "P553", "P554", "P555", "P556", "P557", "P558", + "P559", "P56", "P560", "P561", "P562", "P563", "P564", "P565", + "P566", "P567", "P568", "P569", "P57", "P570", "P571", "P572", + "P573", "P574", "P575", "P576", "P577", "P578", "P579", "P58", + "P580", "P581", "P582", "P583", "P584", "P585", "P586", "P587", + "P588", "P589", "P59", "P590", "P591", "P592", "P593", "P594", + "P595", "P596", "P597", "P598", "P599", "P6", "P60", "P600", + "P601", "P602", "P603", "P604", "P605", "P606", "P607", "P608", + "P609", "P61", "P610", "P611", "P612", "P613", "P614", "P615", + "P616", "P617", "P618", "P619", "P62", "P620", "P621", "P622", + "P623", "P624", "P625", "P626", "P627", "P628", "P629", "P63", + "P630", "P631", "P632", "P633", "P634", "P635", "P636", "P637", + "P638", "P639", "P64", "P640", "P641", "P642", "P643", "P644", + "P645", "P646", "P647", "P648", "P649", "P65", "P650", "P651", + "P652", "P653", "P654", "P655", "P656", "P657", "P658", "P659", + "P66", "P660", "P661", "P662", "P663", "P664", "P665", "P666", + "P667", "P668", "P669", "P67", "P670", "P671", "P672", "P673", + "P674", "P675", "P676", "P677", "P678", "P679", "P68", "P680", + "P681", "P682", "P683", "P684", "P685", "P686", "P687", "P688", + "P689", "P69", "P690", "P691", "P692", "P693", "P694", "P695", + "P696", "P697", "P698", "P699", "P7", "P70", "P700", "P701", + "P702", "P703", "P704", "P705", "P706", "P707", "P708", "P709", + "P71", "P710", "P711", "P712", "P713", "P714", "P715", "P716", + "P717", "P718", "P719", "P72", "P720", "P721", "P722", "P723", + "P724", "P725", "P726", "P727", "P728", "P729", "P73", "P730", + "P731", "P732", "P733", "P734", "P735", "P736", "P737", "P738", + "P739", "P74", "P740", "P741", "P742", "P743", "P744", "P745", + "P746", "P747", "P748", "P749", "P75", "P750", "P751", "P752", + "P753", "P754", "P755", "P756", "P757", "P758", "P759", "P76", + "P760", "P761", "P762", "P763", "P764", "P765", "P766", "P767", + "P768", "P769", "P77", "P770", "P771", "P772", "P773", "P774", + "P775", "P776", "P777", "P778", "P779", "P78", "P780", "P781", + "P782", "P783", "P784", "P785", "P786", "P787", "P788", "P789", + "P79", "P790", "P791", "P792", "P793", "P794", "P795", "P796", + "P797", "P798", "P799", "P8", "P80", "P800", "P801", "P802", + "P803", "P804", "P805", "P806", "P807", "P808", "P809", "P81", + "P810", "P811", "P812", "P813", "P814", "P815", "P816", "P817", + "P818", "P819", "P82", "P820", "P821", "P822", "P823", "P824", + "P825", "P826", "P827", "P828", "P829", "P83", "P830", "P831", + "P832", "P833", "P834", "P835", "P836", "P837", "P838", "P839", + "P84", "P840", "P841", "P842", "P843", "P844", "P845", "P846", + "P847", "P848", "P849", "P85", "P850", "P851", "P852", "P853", + "P854", "P855", "P856", "P857", "P858", "P859", "P86", "P860", + "P861", "P862", "P863", "P864", "P865", "P866", "P867", "P868", + "P869", "P87", "P870", "P871", "P872", "P873", "P874", "P875", + "P876", "P877", "P878", "P879", "P88", "P880", "P881", "P882", + "P883", "P884", "P885", "P886", "P887", "P888", "P889", "P89", + "P890", "P891", "P892", "P893", "P894", "P895", "P896", "P897", + "P898", "P899", "P9", "P90", "P900", "P901", "P902", "P903", + "P904", "P905", "P906", "P907", "P908", "P909", "P91", "P910", + "P911", "P912", "P913", "P914", "P915", "P916", "P917", "P918", + "P919", "P92", "P920", "P921", "P922", "P923", "P924", "P925", + "P926", "P927", "P928", "P929", "P93", "P930", "P931", "P932", + "P933", "P934", "P935", "P936", "P937", "P938", "P939", "P94", + "P940", "P941", "P942", "P943", "P944", "P945", "P946", "P947", + "P948", "P949", "P95", "P950", "P951", "P952", "P953", "P954", + "P955", "P956", "P957", "P958", "P959", "P96", "P960", "P961", + "P962", "P963", "P964", "P965", "P966", "P967", "P968", "P969", + "P97", "P970", "P971", "P972", "P973", "P974", "P975", "P976", + "P977", "P978", "P979", "P98", "P980", "P981", "P982", "P983", + "P984", "P985", "P986", "P987", "P988", "P99") + diff --git a/tests/testthat/_snaps/df_to_array.md b/tests/testthat/_snaps/df_to_array.md new file mode 100644 index 00000000..b2c92f13 --- /dev/null +++ b/tests/testthat/_snaps/df_to_array.md @@ -0,0 +1,111 @@ +# df_to_array() produces consistent results + + WAoAAAACAAQEAQACAwAAAAMOAAACWEAUZmZmZmZmQAwAAAAAAAA/yZmZmZmZmj/2ZmZmZmZm + QBwAAAAAAABACZmZmZmZmj/2ZmZmZmZmQBLMzMzMzM1AGTMzMzMzM0AKZmZmZmZmQAQAAAAA + AABAGAAAAAAAAEATmZmZmZmaQAgAAAAAAAA/yZmZmZmZmj/2ZmZmZmZmQBmZmZmZmZpACZmZ + mZmZmj/4AAAAAAAAQBIAAAAAAABAFzMzMzMzM0AFmZmZmZmaP/5mZmZmZmZAFGZmZmZmZkAS + zMzMzMzNQAmZmZmZmZo/yZmZmZmZmj/0zMzMzMzNQBuZmZmZmZpACMzMzMzMzT/4AAAAAAAA + QBOZmZmZmZpAHGZmZmZmZkAIAAAAAAAAQADMzMzMzM1AF5mZmZmZmkASZmZmZmZmQAjMzMzM + zM0/yZmZmZmZmj/4AAAAAAAAQBYAAAAAAABAAmZmZmZmZj/0zMzMzMzNQBAAAAAAAABAGTMz + MzMzM0AHMzMzMzMzP/zMzMzMzM1AFmZmZmZmZkAUAAAAAAAAQAzMzMzMzM0/yZmZmZmZmj/2 + ZmZmZmZmQBoAAAAAAABABmZmZmZmZj/4AAAAAAAAQBJmZmZmZmZAGgAAAAAAAEAIAAAAAAAA + QAGZmZmZmZpAFzMzMzMzM0AVmZmZmZmaQA8zMzMzMzM/2ZmZmZmZmj/7MzMzMzMzQBbMzMzM + zM1ABmZmZmZmZj/0zMzMzMzNQBIAAAAAAABAHmZmZmZmZkAIAAAAAAAAQADMzMzMzM1AGmZm + ZmZmZkASZmZmZmZmQAszMzMzMzM/0zMzMzMzMz/2ZmZmZmZmQBkzMzMzMzNACmZmZmZmZj/5 + mZmZmZmaQBLMzMzMzM1AE5mZmZmZmkAEAAAAAAAAP/szMzMzMzNAEgAAAAAAAEAUAAAAAAAA + QAszMzMzMzM/yZmZmZmZmj/4AAAAAAAAQBOZmZmZmZpAAzMzMzMzMz/wAAAAAAAAQApmZmZm + ZmZAHTMzMzMzM0AHMzMzMzMzP/zMzMzMzM1AGTMzMzMzM0ARmZmZmZmaQAczMzMzMzM/yZmZ + mZmZmj/2ZmZmZmZmQBpmZmZmZmZABzMzMzMzMz/0zMzMzMzNQBJmZmZmZmZAGszMzMzMzUAE + AAAAAAAAP/zMzMzMzM1AFzMzMzMzM0ATmZmZmZmaQAjMzMzMzM0/uZmZmZmZmj/4AAAAAAAA + QBTMzMzMzM1ABZmZmZmZmj/2ZmZmZmZmQA8zMzMzMzNAHMzMzMzMzUAMzMzMzMzNQAQAAAAA + AABAGGZmZmZmZkAVmZmZmZmaQA2ZmZmZmZo/yZmZmZmZmj/4AAAAAAAAQBQAAAAAAABAAAAA + AAAAAD/wAAAAAAAAQAwAAAAAAABAGgAAAAAAAEAJmZmZmZmaQAAAAAAAAABAFGZmZmZmZkAT + MzMzMzMzQAszMzMzMzM/yZmZmZmZmj/5mZmZmZmaQBeZmZmZmZpACAAAAAAAAD/4AAAAAAAA + QBDMzMzMzM1AGZmZmZmZmkAFmZmZmZmaP/5mZmZmZmZAFTMzMzMzM0ATMzMzMzMzQAgAAAAA + AAA/uZmZmZmZmj/2ZmZmZmZmQBgAAAAAAABAAZmZmZmZmj/wAAAAAAAAQBAAAAAAAABAGzMz + MzMzM0AIAAAAAAAAQADMzMzMzM1AFgAAAAAAAEARMzMzMzMzQAgAAAAAAAA/uZmZmZmZmj/x + mZmZmZmaQBhmZmZmZmZABzMzMzMzMz/2ZmZmZmZmQBLMzMzMzM1AFszMzMzMzUAEAAAAAAAA + QAAAAAAAAABAFAAAAAAAAEAXMzMzMzMzQBAAAAAAAAA/yZmZmZmZmj/zMzMzMzMzQBZmZmZm + ZmZABzMzMzMzMz/0zMzMzMzNQAzMzMzMzM1AFzMzMzMzM0AGZmZmZmZmQAMzMzMzMzNAFGZm + ZmZmZkAWzMzMzMzNQBGZmZmZmZo/2ZmZmZmZmj/4AAAAAAAAQBrMzMzMzM1ACMzMzMzMzT/2 + ZmZmZmZmQBGZmZmZmZpAGZmZmZmZmkAJmZmZmZmaQAJmZmZmZmZAFTMzMzMzM0AVmZmZmZma + QA8zMzMzMzM/2ZmZmZmZmj/0zMzMzMzNQBZmZmZmZmZACAAAAAAAAD/4AAAAAAAAQBIAAAAA + AABAGgAAAAAAAEAIAAAAAAAAP/zMzMzMzM1AFgAAAAAAAEAUZmZmZmZmQAwAAAAAAAA/0zMz + MzMzMz/2ZmZmZmZmQBczMzMzMzNABZmZmZmZmj/wAAAAAAAAQBBmZmZmZmZAHszMzMzMzUAO + ZmZmZmZmQAGZmZmZmZpAGszMzMzMzUAWzMzMzMzNQA5mZmZmZmY/0zMzMzMzMz/7MzMzMzMz + QBjMzMzMzM1AAZmZmZmZmj/4AAAAAAAAQBIAAAAAAABAHszMzMzMzUAEzMzMzMzNQAJmZmZm + ZmZAG5mZmZmZmkAUZmZmZmZmQA5mZmZmZmY/0zMzMzMzMz/4AAAAAAAAQBZmZmZmZmZABAAA + AAAAAD/xmZmZmZmaQA8zMzMzMzNAGAAAAAAAAEABmZmZmZmaP/gAAAAAAABAFAAAAAAAAEAV + mZmZmZmaQAszMzMzMzM/yZmZmZmZmj/7MzMzMzMzQBeZmZmZmZpACZmZmZmZmj/8zMzMzMzN + QBMzMzMzMzNAG5mZmZmZmkAJmZmZmZmaQAJmZmZmZmZAFszMzMzMzUAUZmZmZmZmQA2ZmZmZ + mZo/2ZmZmZmZmj/4AAAAAAAAQBhmZmZmZmZABmZmZmZmZj/0zMzMzMzNQBAAAAAAAABAFmZm + ZmZmZkAGZmZmZmZmQAAAAAAAAABAE5mZmZmZmkASZmZmZmZmQAzMzMzMzM0/yZmZmZmZmj/w + AAAAAAAAQBkzMzMzMzNABAAAAAAAAD/4AAAAAAAAQBOZmZmZmZpAHszMzMzMzUAGZmZmZmZm + QAAAAAAAAABAGszMzMzMzUAUZmZmZmZmQApmZmZmZmY/4AAAAAAAAD/7MzMzMzMzQBhmZmZm + ZmZABmZmZmZmZj/zMzMzMzMzQBLMzMzMzM1AGTMzMzMzM0AFmZmZmZmaP/zMzMzMzM1AE5mZ + mZmZmkATMzMzMzMzQAszMzMzMzM/yZmZmZmZmj/+ZmZmZmZmQBmZmZmZmZpABzMzMzMzMz/0 + zMzMzMzNQBEzMzMzMzNAGszMzMzMzUAKZmZmZmZmQADMzMzMzM1AFszMzMzMzUAUAAAAAAAA + QAgAAAAAAAA/yZmZmZmZmj/5mZmZmZmaQBpmZmZmZmZACAAAAAAAAD/2ZmZmZmZmQBGZmZmZ + mZpAHMzMzMzMzUAJmZmZmZmaP/zMzMzMzM1AGAAAAAAAAEAUAAAAAAAAQAszMzMzMzM/2ZmZ + mZmZmj/5mZmZmZmaQBszMzMzMzNABmZmZmZmZj/2ZmZmZmZmQBMzMzMzMzNAGMzMzMzMzUAG + ZmZmZmZmP/zMzMzMzM1AEzMzMzMzM0AUzMzMzMzNQAwAAAAAAAA/yZmZmZmZmj/4AAAAAAAA + QBrMzMzMzM1ACAAAAAAAAD/7MzMzMzMzQBQAAAAAAABAGGZmZmZmZkAIAAAAAAAAP/zMzMzM + zM1AE5mZmZmZmkAUzMzMzMzNQAszMzMzMzM/yZmZmZmZmj/2ZmZmZmZmQBgAAAAAAABABzMz + MzMzMz/4AAAAAAAAQBIAAAAAAABAGZmZmZmZmkAGZmZmZmZmQADMzMzMzM1AFmZmZmZmZkAS + zMzMzMzNQAmZmZmZmZo/yZmZmZmZmj/5mZmZmZmaQBbMzMzMzM1ABMzMzMzMzT/wAAAAAAAA + QAwAAAAAAABAHMzMzMzMzUAIAAAAAAAAP/mZmZmZmZpAFzMzMzMzM0ATMzMzMzMzQAjMzMzM + zM0/yZmZmZmZmj/5mZmZmZmaQBYAAAAAAABAAzMzMzMzMz/xmZmZmZmaQA5mZmZmZmZAHZmZ + mZmZmkAGZmZmZmZmP/5mZmZmZmZAGGZmZmZmZkAVmZmZmZmaQAszMzMzMzM/2ZmZmZmZmj/4 + AAAAAAAAQBYAAAAAAABAAzMzMzMzMz/wAAAAAAAAQA2ZmZmZmZpAH5mZmZmZmkAOZmZmZmZm + QAAAAAAAAABAGZmZmZmZmkAUzMzMzMzNQBBmZmZmZmY/uZmZmZmZmj/4AAAAAAAAQBczMzMz + MzNABZmZmZmZmj/zMzMzMzMzQA8zMzMzMzNAGZmZmZmZmkAGZmZmZmZmQAGZmZmZmZpAFmZm + ZmZmZkAWAAAAAAAAQBDMzMzMzM0/yZmZmZmZmj/2ZmZmZmZmQBgAAAAAAABABZmZmZmZmj/5 + mZmZmZmaQBRmZmZmZmZAGTMzMzMzM0AGZmZmZmZmP/gAAAAAAABAFGZmZmZmZkATmZmZmZma + QAjMzMzMzM0/yZmZmZmZmj/4AAAAAAAAQBWZmZmZmZpACAAAAAAAAD/4AAAAAAAAQBIAAAAA + AABAGGZmZmZmZkAEzMzMzMzNP/ZmZmZmZmZAFmZmZmZmZkAUAAAAAAAAQAmZmZmZmZo/yZmZ + mZmZmj/zMzMzMzMzQBgAAAAAAABACzMzMzMzMz/5mZmZmZmaQBIAAAAAAABAHszMzMzMzUAI + AAAAAAAAQAJmZmZmZmZAGGZmZmZmZkAWAAAAAAAAQAwAAAAAAAA/yZmZmZmZmj/0zMzMzMzN + QBrMzMzMzM1ACMzMzMzMzT/4AAAAAAAAQBLMzMzMzM1AGTMzMzMzM0ALMzMzMzMzQAMzMzMz + MzNAFmZmZmZmZkATmZmZmZmaQAzMzMzMzM0/uZmZmZmZmj/2ZmZmZmZmQBkzMzMzMzNAAmZm + ZmZmZj/0zMzMzMzNQBGZmZmZmZpAGZmZmZmZmkAIzMzMzMzNP/zMzMzMzM1AFgAAAAAAAEAR + mZmZmZmaQAgAAAAAAAA/yZmZmZmZmj/0zMzMzMzNQBZmZmZmZmZACAAAAAAAAD/0zMzMzMzN + QBBmZmZmZmZAGAAAAAAAAEAIAAAAAAAAP/zMzMzMzM1AEzMzMzMzM0AUZmZmZmZmQAszMzMz + MzM/yZmZmZmZmj/4AAAAAAAAQBYAAAAAAABABAAAAAAAAD/0zMzMzMzNQBAAAAAAAABAG5mZ + mZmZmkAIzMzMzMzNQADMzMzMzM1AFZmZmZmZmkAUAAAAAAAAQAwAAAAAAAA/0zMzMzMzMz/0 + zMzMzMzNQBYAAAAAAABABMzMzMzMzT/zMzMzMzMzQBGZmZmZmZpAGszMzMzMzUAIzMzMzMzN + QAMzMzMzMzNAFmZmZmZmZkASAAAAAAAAQAJmZmZmZmY/0zMzMzMzMz/0zMzMzMzNQBhmZmZm + ZmZACAAAAAAAAD/2ZmZmZmZmQBJmZmZmZmZAG5mZmZmZmkAIzMzMzMzNQAJmZmZmZmZAFGZm + ZmZmZkARmZmZmZmaQAmZmZmZmZo/yZmZmZmZmj/0zMzMzMzNQBczMzMzMzNABMzMzMzMzT/z + MzMzMzMzQBAAAAAAAABAFzMzMzMzM0AFmZmZmZmaP/5mZmZmZmZAFGZmZmZmZkAUAAAAAAAA + QAwAAAAAAAA/4zMzMzMzMz/5mZmZmZmaQBQAAAAAAABAAmZmZmZmZj/wAAAAAAAAQApmZmZm + ZmZAGzMzMzMzM0AJmZmZmZmaQAJmZmZmZmZAF5mZmZmZmkAUZmZmZmZmQA5mZmZmZmY/2ZmZ + mZmZmj/+ZmZmZmZmQBZmZmZmZmZABZmZmZmZmj/0zMzMzMzNQBDMzMzMzM1AGszMzMzMzUAK + ZmZmZmZmQAQAAAAAAABAFszMzMzMzUATMzMzMzMzQAgAAAAAAAA/0zMzMzMzMz/2ZmZmZmZm + QBbMzMzMzM1ACAAAAAAAAD/zMzMzMzMzQBDMzMzMzM1AGszMzMzMzUAIAAAAAAAAQAJmZmZm + ZmZAFMzMzMzMzUAUZmZmZmZmQA5mZmZmZmY/yZmZmZmZmj/5mZmZmZmaQBbMzMzMzM1ABzMz + MzMzMz/0zMzMzMzNQBDMzMzMzM1AGTMzMzMzM0AEAAAAAAAAP/5mZmZmZmZAFAAAAAAAAEAS + ZmZmZmZmQAmZmZmZmZo/yZmZmZmZmj/2ZmZmZmZmQBjMzMzMzM1ABzMzMzMzMz/0zMzMzMzN + QBEzMzMzMzNAGgAAAAAAAEAIAAAAAAAAQAAAAAAAAABAFMzMzMzMzUAVMzMzMzMzQA2ZmZmZ + mZo/yZmZmZmZmj/4AAAAAAAAQBRmZmZmZmZABAAAAAAAAD/xmZmZmZmaQAgAAAAAAABAGMzM + zMzMzUALMzMzMzMzQAJmZmZmZmZAFZmZmZmZmkAUAAAAAAAAQApmZmZmZmY/yZmZmZmZmj/2 + ZmZmZmZmQBbMzMzMzM1ABmZmZmZmZj/0zMzMzMzNQBBmZmZmZmZAF5mZmZmZmkAIAAAAAAAA + P/zMzMzMzM1AFGZmZmZmZgAABAIAAAABAAQACQAAAANkaW0AAAANAAAAAwAAAAQAAAADAAAA + MgAABAIAAAABAAQACQAAAAhkaW1uYW1lcwAAAhMAAAADAAAAEAAAAAQABAAJAAAADFNlcGFs + Lkxlbmd0aAAEAAkAAAALU2VwYWwuV2lkdGgABAAJAAAAC1BldGFsLldpZHRoAAQACQAAAAxQ + ZXRhbC5MZW5ndGgAAAAQAAAAAwAEAAkAAAAGc2V0b3NhAAQACQAAAAp2ZXJzaWNvbG9yAAQA + CQAAAAl2aXJnaW5pY2EAAAAQAAAAMgAEAAkAAAABMQAEAAkAAAABMgAEAAkAAAABMwAEAAkA + AAABNAAEAAkAAAABNQAEAAkAAAABNgAEAAkAAAABNwAEAAkAAAABOAAEAAkAAAABOQAEAAkA + AAACMTAABAAJAAAAAjExAAQACQAAAAIxMgAEAAkAAAACMTMABAAJAAAAAjE0AAQACQAAAAIx + NQAEAAkAAAACMTYABAAJAAAAAjE3AAQACQAAAAIxOAAEAAkAAAACMTkABAAJAAAAAjIwAAQA + CQAAAAIyMQAEAAkAAAACMjIABAAJAAAAAjIzAAQACQAAAAIyNAAEAAkAAAACMjUABAAJAAAA + AjI2AAQACQAAAAIyNwAEAAkAAAACMjgABAAJAAAAAjI5AAQACQAAAAIzMAAEAAkAAAACMzEA + BAAJAAAAAjMyAAQACQAAAAIzMwAEAAkAAAACMzQABAAJAAAAAjM1AAQACQAAAAIzNgAEAAkA + AAACMzcABAAJAAAAAjM4AAQACQAAAAIzOQAEAAkAAAACNDAABAAJAAAAAjQxAAQACQAAAAI0 + MgAEAAkAAAACNDMABAAJAAAAAjQ0AAQACQAAAAI0NQAEAAkAAAACNDYABAAJAAAAAjQ3AAQA + CQAAAAI0OAAEAAkAAAACNDkABAAJAAAAAjUwAAAEAgAAAAEABAAJAAAABW5hbWVzAAAAEAAA + AAMABAAJAAAACXBhcmFtZXRlcgAEAAkAAAAHU3BlY2llcwAEAAkAAAADb2JzAAAA/gAABAIA + AAABAAQACQAAAAVjbGFzcwAAABAAAAACAAQACQAAAAV4dGFicwAEAAkAAAAFdGFibGUAAAQC + AAAAAQAEAAkAAAAEY2FsbAAAAAYAAAABAAQACQAAAAV4dGFicwAABAIAAAABAAQACQAAAAdm + b3JtdWxhAAAABgAAB/8AAAACAAAAAQAEAAkAAAANeHRhYnNfZm9ybXVsYQAAAP4AAAQCAAAA + AQAEAAkAAAAEZGF0YQAAAAEABAAJAAAAAS4AAAD+AAAA/g== + diff --git a/tests/testthat/_snaps/log_likelihood.md b/tests/testthat/_snaps/log_likelihood.md new file mode 100644 index 00000000..4fd146e4 --- /dev/null +++ b/tests/testthat/_snaps/log_likelihood.md @@ -0,0 +1,4 @@ +# `log_likelihood()` gives consistent results + + -9268.8238 + diff --git a/tests/testthat/test-as_curve_params.R b/tests/testthat/test-as_curve_params.R index c89f457d..6e5e6a83 100644 --- a/tests/testthat/test-as_curve_params.R +++ b/tests/testthat/test-as_curve_params.R @@ -1,15 +1,17 @@ -test_that("`as_curve_params()` produces an error when non-curve data is provided", { - library(magrittr) - expect_error( - object = curve_data <- - "https://osf.io/download//n6cp3/" %>% # pop data - readr::read_rds() %>% - as_curve_params(), - class = "not curve_params" - ) -}) +test_that("`as_curve_params()` produces an error + when non-curve data is provided", { + library(magrittr) + expect_error( + object = curve_data <- + "https://osf.io/download//n6cp3/" %>% # pop data + readr::read_rds() %>% + as_curve_params(), + class = "not curve_params" + ) + }) -test_that("`as_curve_params()` produces an error when `data` is not a data.frame", +test_that("`as_curve_params()` produces an error + when `data` is not a data.frame", { library(magrittr) expect_error(object = @@ -28,5 +30,7 @@ test_that("`as_curve_params()` produces expected results", { expect_snapshot_value(x = test_data, style = "serialize") + test_data %>% ssdtools:::expect_snapshot_data(name = "curve-data") + }) diff --git a/tests/testthat/test-class_attributes.R b/tests/testthat/test-class_attributes.R new file mode 100644 index 00000000..a94efd7b --- /dev/null +++ b/tests/testthat/test-class_attributes.R @@ -0,0 +1,36 @@ +test_that("`get_biomarker_levels()` works", { + xs_data <- "https://osf.io/download//n6cp3/" %>% + load_pop_data() + biomarker_levels <- xs_data %>% get_biomarker_levels() + expected_levels <- structure(1:2, + levels = c("HlyE_IgA", "HlyE_IgG"), + class = "factor") + expect_equal(object = biomarker_levels, expected = expected_levels) +}) + +test_that("`get_id()` works", { + xs_data <- "https://osf.io/download//n6cp3/" %>% + load_pop_data() + + xs_data %>% + get_id() %>% + sort() %>% + expect_snapshot_value(style = "deparse") + +}) + +test_that("`get_biomarker_names_var() works", { + biomarker_names_var <- + "https://osf.io/download//n6cp3/" %>% + load_pop_data() %>% + get_biomarker_names_var() + + expect_equal(object = biomarker_names_var, expected = "antigen_iso") +}) + + +test_that("`set_age()` detects partial matches", { + "https://osf.io/download//n6cp3/" %>% + load_pop_data(age = "age$") %>% + expect_warning(class = "missing variable") +}) diff --git a/tests/testthat/test-df_to_array.R b/tests/testthat/test-df_to_array.R new file mode 100644 index 00000000..39b760aa --- /dev/null +++ b/tests/testthat/test-df_to_array.R @@ -0,0 +1,18 @@ +test_that("df_to_array() produces consistent results", { + library(dplyr) + library(tidyr) + df <- iris %>% + tidyr::pivot_longer( + names_to = "parameter", + cols = c( + "Sepal.Length", + "Sepal.Width", + "Petal.Width", + "Petal.Length" + ) + ) %>% + mutate(parameter = factor(parameter, levels = unique(parameter))) + arr <- df %>% + serocalculator:::df_to_array(dim_var_names = c("parameter", "Species")) + arr %>% expect_snapshot_value(style = "serialize") +}) diff --git a/tests/testthat/test-log_likelihood.R b/tests/testthat/test-log_likelihood.R new file mode 100644 index 00000000..0a4e3d72 --- /dev/null +++ b/tests/testthat/test-log_likelihood.R @@ -0,0 +1,34 @@ +test_that("`log_likelihood()` gives consistent results", { + library(dplyr) + library(tibble) + + # load in longitudinal parameters + dmcmc <- load_curve_params("https://osf.io/download/rtw5k") + + xs_data <- "https://osf.io/download//n6cp3/" %>% + load_pop_data() + + # Load noise params + cond <- tibble( + antigen_iso = c("HlyE_IgG", "HlyE_IgA"), + nu = c(0.5, 0.5), + # Biologic noise (nu) + eps = c(0, 0), + # M noise (eps) + y.low = c(1, 1), + # low cutoff (llod) + y.high = c(5e6, 5e6) + ) # high cutoff (y.high) + + # Calculate log-likelihood + ll_AG <- log_likelihood( # nolint: object_name_linter + pop_data = xs_data, + curve_params = dmcmc, + noise_params = cond, + antigen_isos = c("HlyE_IgG", "HlyE_IgA"), + lambda = 0.1 + ) + + expect_snapshot_value(ll_AG) + +}) diff --git a/tests/testthat/test-set_biomarker_var.R b/tests/testthat/test-set_biomarker_var.R new file mode 100644 index 00000000..31feb13b --- /dev/null +++ b/tests/testthat/test-set_biomarker_var.R @@ -0,0 +1,6 @@ +test_that("`set_biomarker_var()` halts on misspecified column", { + xs_data <- "https://osf.io/download//n6cp3/" %>% + readr::read_rds() %>% + set_biomarker_var("biomarker") %>% + expect_error(class = "missing variable") +}) diff --git a/vignettes/.gitignore b/vignettes/.gitignore index 9289a5e3..05becfb6 100644 --- a/vignettes/.gitignore +++ b/vignettes/.gitignore @@ -2,3 +2,10 @@ *.R methodology.log methodology_files +~$*.pptx +~$*.docx +*_files/ +*.pptx +*.docx + +/.quarto/ diff --git a/vignettes/articles/_antibody-response-model.qmd b/vignettes/articles/_antibody-response-model.qmd new file mode 100644 index 00000000..d867a851 --- /dev/null +++ b/vignettes/articles/_antibody-response-model.qmd @@ -0,0 +1,221 @@ +# Modeling the seroresponse kinetics curve + +::: notes +Now, we need a model for the antibody response to infection, $\pdf(Y=y|T=t)$. +The current version of the `{serocalculator}` package uses +a two-phase model for +the shape of the seroresponse [@Teunis_2016]. +::: + +## Model for active infection period {.smaller} + +::: notes +The first phase of the model represents the **active infection period**, +and uses a simplified Lotka-Volterra predator-prey model +[@volterra1928variations] +where the pathogen is the prey and the antibodies are the predator: +::: + +Notation: + +* $x(t)$: Pathogen concentration at time $t$ +* $y(t)$: Antibody concentration at time $t$ + +Model: + +* $$x'(t) = \alpha x(t) - \beta y(t)$$ +* $$y'(t) = \delta y(t)$$ + +::: notes +With baseline antibody concentration $y(0) = y_{0}$ and +initial pathogen concentration +$x(0) = x_{0}$. +::: + +::: notes +Compared to the standard LV model: + +* the predation term with the $\beta$ coefficient is missing +the prey concentration $x(t)$ factor; +we assume that the efficiency of predation +doesn't depend on pathogen concentration. + +* the differential equation for predator density is missing +the predator death rate term +$-\gamma y(t)$; we assume that as long as there are any pathogens present, +the antibody decay rate is negligible compared to the growth rate. + +* the predator growth rate term $\delta y(t)$ is missing +the prey density factor $x(t)$ +we assume that as long as there are any pathogens present, +the antibody concentration grows at the same exponential rate. + +These omissions were made to simplify the estimation process, +under the assumption that they are negligible +compared to the other terms in the model. + +::: + +## Model for post-infection antibody decay + +* $$b(t) = 0$$ + +* $$y^{\prime}(t) = -\alpha y(t)^r$$ + +::: notes +Antibody decay is different from exponential (log--linear) decay. +When the shape parameter $r > 1$, +log concentrations decrease rapidly after infection has terminated, +but decay then slows down and +low antibody concentrations are maintained for a long period. +If $r = 1$, +this model reduces to exponential decay with decay rate $\alpha$. +::: + +## Putting it all together + +The serum antibody response $y(t)$ can be written as + +$$ +y(t) = y_{+}(t) + y_{-}(t) +$$ + +where + +$$ +\begin{align} +y_{+}(t) & = y_{0}\text{e}^{\mu t}[0\leq t % + load_curve_params() %>% + filter(iter < 50) + +curve1 = + curves %>% + filter( + # iter %in% 1:10, + iter == 5, + antigen_iso == cur_ai) + +library(ggplot2) + +curve1 %>% +serocalculator:::plot_curve_params_one_ab( + log_y = FALSE +) + + xlim(0, 100) + + theme_minimal() + + geom_vline( + aes(xintercept = curve1$t1, + col = "t1") + ) + + + geom_hline( + aes(yintercept = curve1$y0, + col = "y0") + ) + + + + geom_hline( + aes(yintercept = curve1$y1, + col = "y1") + ) + + # geom_point( + # data = curve1, + # aes( + # x = t1, + # y = y1, + # col = "(t1,y1)" + # ) + # ) + + theme(legend.position = "bottom") + + labs(col = "") + +``` + + +:::: notes +The antibody level at $t=0$ is $y_{0}$; +the rising branch ends at $t = t_{1}$ +where the peak antibody level $y_{1}$ is reached. +Any antibody level $y(t) \in (y_{0}, y_{1})$ eventually occurs twice. +:::: + +--- + +::: notes +An interactive Shiny app that allows the user to manipulate the model parameters is available at: + + + +::: + +```{r} +#| echo: false +#| eval: false +#| label: qr-code +library(qrcode) +code <- qr_code("https://ucdserg.shinyapps.io/antibody-kinetics-model-2/") +generate_svg(code, filename = "man/figures/qr_akm2.svg") +``` + +![QR code for shiny app](../man/figures/qr_akm2.svg) + +## Biological noise + +::: notes +When we measure antibody concentrations in a blood sample, we are essentially counting molecules (using biochemistry). + +We might miss some of the antibodies (undercount, false negatives) and we also might incorrectly count some other molecules that aren't actually the ones we are looking for (overcount, false positives, cross-reactivity). + +We are more concerned about overcount (cross-reactivity) than undercount. For a given antibody, we can do some analytical work beforehand to estimate the distribution of overcounts, and add that to our model $p(Y=y|T=t)$. +::: + +Notation: + +* $y_\text{obs}$: measured serum antibody concentration +* $y_\text{true}$: "true" serum antibody concentration +* $\epsilon_b$: noise due to probe cross-reactivity + +Model: + +* $y_\text{obs} = y_\text{true} + \epsilon_b$ +* $\epsilon_b \sim \text{Unif}(0, \nu)$ + +$\nu$ needs to be pre-estimated using negative controls, +typically using the 95th percentile of the distribution of +antibody responses to the antigen-isotype +in a population with no exposure. + +## Measurement noise + +There are also some other sources of noise in our bioassays; +user differences in pipetting technique, random ELISA plate effects, etc. +This noise can cause both overcount and undercount. We can also estimate the magnitude of this noise source, and include it in $p(Y=y|T=t)$. + +Measurement noise, $\varepsilon$ ("epsilon"), represents measurement error from the laboratory testing process. +It is defined by a CV (coefficient of variation) as the ratio of the standard deviation to the mean for replicates. +Note that the CV should ideally be measured across plates rather than within the same plate. diff --git a/vignettes/articles/_macros.qmd b/vignettes/articles/_macros.qmd new file mode 100644 index 00000000..269684e0 --- /dev/null +++ b/vignettes/articles/_macros.qmd @@ -0,0 +1,9 @@ +\def\NA{\text{NA}} +\def\pdf{\mathbb{p}} +\def\pmf{\mathbb{P}} +\def\Lik{\mathcal{L}} +\def\llik{\ell} +\renewcommand{\vec}[1]{\mathbf{#1}} +\providecommand{\expf}[1]{\exp{\left\{#1\right\}}} +\providecommand{\logf}[1]{\log{\left\{#1\right\}}} +\providecommand{\invf}[1]{#1^{-1}} diff --git a/vignettes/articles/_methods-continued.qmd b/vignettes/articles/_methods-continued.qmd new file mode 100644 index 00000000..4887463f --- /dev/null +++ b/vignettes/articles/_methods-continued.qmd @@ -0,0 +1,291 @@ + +# Other innovations in the `serocalculator` package + +The `serocalculator` package provides three refinements to the method for +calculating seroincidence published earlier [@Teunis_2012] and implemented in R package +`seroincidence`, versions 1.x: + +(1) inclusion of infection episode with rising antibody levels, + +(2) non--exponential decay of serum antibodies after infection, and + +(3) age--dependent correction for subjects who have never seroconverted. + +It is important to note that, although the implemented +methods use a specific parametric model, as proposed in [@de_Graaf_2014] and augmented in +[@Teunis_2016], the methods used to calculate the likelihood function allow +seroresponses of arbitrary shape. + + +# Backward recurrence time + +Considering the (fundamental) uniform distribution $u_{f}$ of sampling within a given interval, +the interval length distribution $p(\Delta t)$ +and the distribution of (cross-sectionally) sampled +interval length [@Teunis_2012] + +$$ +q(\Delta t) = \frac{p(\Delta t)\cdot\Delta t}{\overline{\Delta t_{p}}} +$$ +the joint distribution of backward recurrence time +and cross--sectional interval length +is the product $u_{f}\cdot q$, +because these probabilities are independent. + +The distribution of backward recurrence time is the marginal distribution + +$$ +\begin{align} +u(\tau) & = \int_{\Delta t=0}^{\infty} u_{f}(\tau;\Delta t)\cdot q(\Delta t)\text{d}\Delta t\\ + & = \int_{0}^{\infty}\frac{[0\le\tau\le\Delta t]}{\Delta t}\cdot \frac{p(\Delta t)\cdot + \Delta t}{\overline{\Delta t_{p}}}\text{d}\Delta t\\ + & = \frac{1}{\overline{\Delta t_{p}}}\int_{\tau}^{\infty}p(\Delta t)\text{d}\Delta t +\end{align} +$$ + +# Incidence of seroconversions + +To calculate the incidence of seroconversions, as in [@Teunis_2012], the distribution +$p(\Delta t)$ of intervals $\Delta t$ between seroconversions, is important. Assuming any subject is +sampled completely at random during their intervals between seroconversions, and accounting for +interval length bias [@satten_etal04; @zelen04], the distribution of backward recurrence times +$\tau$ can be written as [@Teunis_2012] + +$$ +u(\tau) = \frac{1}{\overline{\Delta t}} \int_{\tau=0}^{\infty}p(\Delta t)\text{d}\Delta t + = \frac{1-P(\Delta t)}{\overline{\Delta t}} +$$ +where $\overline{\Delta t}$ is the $p$--distribution expected value of $\Delta t$. + +This density is employed to estimate seroconversion rates. The antibody concentration $y$ is the +observable quantity, and we need to express the probability (density) of observed $y$ in terms of +the density of backward recurrence time. + +First, the backward recurrence time can $\tau$ be expressed as a function of the serum antibody +concentration $y$ +$$ +\tau(y) = \tau_{+}(y) + \tau_{-}(y) +$$ +where + +$$ +\begin{align} +\tau_{+}(y) & = \displaystyle{\frac{1}{\mu}} \log\left(\displaystyle{\frac{y_{+}}{y_{0}}}\right) [0\le \tau %\VignetteIndexEntry{Enteric Fever Seroincidence Vignette} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} -bibliography: references.bib - --- ## Introduction @@ -161,7 +159,7 @@ Once log transformed, our data looks much more normally distributed. In most cas Let's also take a look at how antibody responses change by age. -```{r "plot age", message=FALSE, warning = FALSE} +```{r plot-age, message=FALSE, warning = FALSE} #Plot antibody responses by age ggplot(data = xs_data %>% diff --git a/vignettes/articles/fig/fig1a-1.svg b/vignettes/articles/fig/fig1a-1.svg new file mode 100644 index 00000000..f68e997e --- /dev/null +++ b/vignettes/articles/fig/fig1a-1.svg @@ -0,0 +1,3795 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/vignettes/articles/fig/response.pdf b/vignettes/articles/fig/response.pdf new file mode 100644 index 00000000..f4a78eb2 Binary files /dev/null and b/vignettes/articles/fig/response.pdf differ diff --git a/vignettes/articles/references.bib b/vignettes/articles/references.bib deleted file mode 100644 index 7b76651a..00000000 --- a/vignettes/articles/references.bib +++ /dev/null @@ -1,109 +0,0 @@ -@article{de_Graaf_2014, - doi = {10.1016/j.epidem.2014.08.002}, - url = {https://doi.org/10.1016%2Fj.epidem.2014.08.002}, - year = 2014, - month = {dec}, - publisher = {Elsevier {BV}}, - volume = {9}, - pages = {1--7}, - author = {WF deGraaf and MEE Kretzschmar and PFM Teunis and O Diekmann}, - title = {A two-phase within-host model for immune response and its application to serological profiles of pertussis}, - journal = {Epidemics} -} - -@article{Simonsen_2009, - doi = {10.1002/sim.3592}, - url = {https://doi.org/10.1002%2Fsim.3592}, - year = 2009, - month = {jun}, - publisher = {Wiley}, - volume = {28}, - number = {14}, - pages = {1882--1895}, - author = {J Simonsen and K M{\o}lbak and G Falkenhorst and KA Krogfelt and A Linneberg and PFM Teunis}, - title = {Estimation of incidences of infectious diseases based on antibody measurements}, - journal = {Statistics in Medicine} -} - -@article{Teunis_2012, - doi = {10.1002/sim.5322}, - url = {https://doi.org/10.1002%2Fsim.5322}, - year = 2012, - month = {mar}, - publisher = {Wiley}, - volume = {31}, - number = {20}, - pages = {2240--2248}, - author = {PFM Teunis and JCH van Eijkeren and CW Ang and YTHP van Duynhoven and JB Simonsen and MA Strid and W van Pelt}, - title = {Biomarker dynamics: estimating infection rates from serological data}, - journal = {Statistics in Medicine} -} - -@article{Teunis_2016, - doi = {10.1016/j.epidem.2016.04.001}, - url = {https://doi.org/10.1016%2Fj.epidem.2016.04.001}, - year = 2016, - month = {sep}, - publisher = {Elsevier {BV}}, - volume = {16}, - pages = {33--39}, - author = {PFM Teunis and JCH van Eijkeren and WF deGraaf and A Bona{\v{c}}i{\'{c}} Marinovi{\'{c}} and MEE Kretzschmar}, - title = {Linking the seroresponse to infection to within-host heterogeneity in antibody production}, - journal = {Epidemics} -} - -@article{Strid_2001, - doi = {10.1128/cdli.8.2.314-319.2001}, - url = {https://doi.org/10.1128%2Fcdli.8.2.314-319.2001}, - year = 2001, - month = {mar}, - publisher = {American Society for Microbiology}, - volume = {8}, - number = {2}, - pages = {314--319}, - author = {MA Strid and J Engberg and LB Larsen and K Begtrup and K M{\o}lbak and KA Krogfelt}, - title = {Antibody responses to Campylobacter infections determined by an enzyme-linked immunosorbent assay: 2-year follow-up study of 210 patients}, - journal = {Clinical Diagnostic Laboratory Immunology} -} - -@article{Aiemjoy_2022_Lancet, - url = {https://doi.org/10.1016/S2666-5247(22)00114-8}, - year = 2022, - volume = {3}, - number = {8}, - pages = {e578--e587}, - author = {K Aiemjoy and JC Seidman and S Saha and SJ Munira and MS Islam Sajib and SM al Sium and A Sarkar and N Alam and FN Zahan and MS Kabir and D Tamrakar and K Vaidya and R Shrestha and J Shakya and N Katuwal and S Shrestha and MT Yousafzai and J Iqbal and IF Dehraj and … and JR Andrews }, - title = {Estimating typhoid incidence from community-based serosurveys: a multicohort study}, - journal = {The Lancet Microbe} -} - -@article{Aiemjoy_2024_scrub, - title={Estimating the seroincidence of scrub typhus using antibody dynamics following infection}, - author={K Aiemjoy and N Katuwal and K Vaidya and S Shrestha and M Thapa and PFM Teunis and II Bogoch and P Trowbridge and P Kantipong and SD Blacksell and others}, - journal={American Journal of Tropical Medicine and Hygiene}, - year={Accepted Feb 2024}, - -} - -@article{Aiemjoy_2022_SouthSudan, - doi = {10.3201/eid2811.220239}, - year = 2022, - month = {nov}, - volume = {28}, - number = {11}, - pages = {2316--2320}, - author = {K Aiemjoy and J Rumunu and JJ Hassen and KE Wiens and D Garrett and P Kamenskaya and JB Harris and AS Azman and PFM Teunis and JC Seidman and JF Wamala and JR Andrews and RC Charles}, - title = {Seroincidence of Enteric Fever,Juba, South Sudan}, - journal = {Emerging Infectious Diseases} -} - -@article{Teunis_2020, - doi = {10.1002/sim.8578}, - year = 2020, - volume = {39}, - number = {21}, - pages = {2799--2814}, - author = {PFM Teunis and JCH van Eijkeren}, - title = {Estimation of seroconversion rates for infectious diseases: Effects of age and noise}, - journal = {Statistics in Medicine} -} diff --git a/vignettes/articles/scrubTyphus_example.Rmd b/vignettes/articles/scrubTyphus_example.Rmd index b1cc5045..1121b518 100644 --- a/vignettes/articles/scrubTyphus_example.Rmd +++ b/vignettes/articles/scrubTyphus_example.Rmd @@ -9,7 +9,6 @@ vignette: > %\VignetteIndexEntry{Scrub Typhus Seroincidence Vignette} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} -bibliography: references.bib --- @@ -121,7 +120,7 @@ xs_data %>% summary(strata = "country") Let's also take a look at how antibody responses change by age. -```{r "plot age"} +```{r plot-age} # Plot antibody responses by age autoplot(object = xs_data, type = "age-scatter", strata = "Country") ``` diff --git a/vignettes/articles/serocalculator.Rmd b/vignettes/articles/serocalculator.Rmd index 76d5ce66..23113bc6 100644 --- a/vignettes/articles/serocalculator.Rmd +++ b/vignettes/articles/serocalculator.Rmd @@ -1,13 +1,12 @@ --- -title: "Methodology" +title: "Overview" description: > A summary of the methods behind serocalculator. output: rmarkdown::html_vignette vignette: > - %\VignetteIndexEntry{Methodology} + %\VignetteIndexEntry{Overview} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} -bibliography: references.bib --- ## Overview diff --git a/vignettes/articles/simulate_xsectionalData.Rmd b/vignettes/articles/simulate_xsectionalData.Rmd index 29273b0c..3869428d 100644 --- a/vignettes/articles/simulate_xsectionalData.Rmd +++ b/vignettes/articles/simulate_xsectionalData.Rmd @@ -11,12 +11,14 @@ vignette: > %\VignetteIndexEntry{simulate_xsectionalData} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} ---- - - + --- + + This vignette shows how to simulate a cross-sectional sample of -seroresponses for incident infections as a Poisson process with frequency `lambda`. -Responses are generated for the antibodies given in the `antigen_isos` argument. +seroresponses for incident infections as a Poisson process with +frequency `lambda`. +Responses are generated for the antibodies given in the `antigen_isos` +argument. Age range of the simulated cross-sectional record is `lifespan`. @@ -30,19 +32,22 @@ However, when `age.fx` is set to NA then the age at infection is used. The boolean `renew.params` determines whether each infection uses a new set of longitudinal parameters, sampled at random from the -posterior predictive output of the longitudinal model. If set to `FALSE` +posterior predictive output of the longitudinal model. +If set to `FALSE`, a parameter set is chosen at birth and kept, but: -1. the baseline antibody levels (`y0`) are updated with the simulated level -(just) prior to infection, and -2. when `is.na(age.fx)` then the selected parameter sample is updated for the -age when infection occurs. +1. the baseline antibody levels (`y0`) are updated with the simulated +level (just) prior to infection, and + +2. when `is.na(age.fx)` then the selected parameter sample is updated +for the age when infection occurs. + There is also a variable `n.mc`: when `n.mc==0` then a random MC sample is chosen out of the posterior set (1:4000). When `n.mc` is given a value in 1:4000 then the chosen number is fixed and reused in any subsequent infection. This is for diagnostic purposes. - - + + ```{r, include = FALSE} knitr::opts_chunk$set( collapse = TRUE, @@ -56,21 +61,25 @@ knitr::opts_chunk$set( ## load model parameters +Here we load in longitudinal parameters; +these are modeled from all SEES cases across all ages and countries: + ```{r setup} library(serocalculator) library(tidyverse) library(ggbeeswarm) # for plotting - -# load in longitudinal parameters, these are modeled from all SEES cases across all ages and countries +library(dplyr) dmcmc <- "https://osf.io/download/rtw5k" %>% load_curve_params() %>% - filter(iter < 500) # reduce number of mcmc samples for speed + dplyr::filter(iter < 500) # reduce number of mcmc samples for speed ``` ## visualize antibody decay model -We can graph individual MCMC samples from the posterior distribution of model parameters using a `autoplot.curve_params()` method for the `autoplot()` function: +We can graph individual MCMC samples from the posterior distribution +of model parameters using a `autoplot.curve_params()` method for the +`autoplot()` function: ```{r} dmcmc %>% autoplot(n_curves = 50) @@ -155,7 +164,12 @@ ggplot(csdata, aes( x = as.factor(antigen_iso), y = value )) + - geom_beeswarm(size = .2, alpha = .3, aes(color = antigen_iso), show.legend = F) + + geom_beeswarm( + size = .2, + alpha = .3, + aes(color = antigen_iso), + show.legend = FALSE + ) + geom_boxplot(outlier.colour = NA, fill = NA) + scale_y_log10() + theme_linedraw() + @@ -167,47 +181,54 @@ ggplot(csdata, aes( We can calculate the log-likelihood of the data as a function of the incidence rate directly: ```{r} -ll_A <- log_likelihood( - pop_data = csdata, - curve_params = dmcmc, - noise_params = cond, - antigen_isos = "HlyE_IgA", - lambda = 0.1 -) %>% print() +ll_a <- + log_likelihood( + pop_data = csdata, + curve_params = dmcmc, + noise_params = cond, + antigen_isos = "HlyE_IgA", + lambda = 0.1 + ) %>% + print() -ll_G <- log_likelihood( - pop_data = csdata, - curve_params = dmcmc, - noise_params = cond, - antigen_isos = "HlyE_IgG", - lambda = 0.1 -) %>% print() +ll_g <- + log_likelihood( + pop_data = csdata, + curve_params = dmcmc, + noise_params = cond, + antigen_isos = "HlyE_IgG", + lambda = 0.1 + ) %>% + print() -ll_AG <- llik( - pop_data = csdata, - curve_params = dmcmc, - noise_params = cond, - antigen_isos = c("HlyE_IgG", "HlyE_IgA"), - lambda = 0.1 -) %>% print() +ll_ag <- + log_likelihood( + pop_data = csdata, + curve_params = dmcmc, + noise_params = cond, + antigen_isos = c("HlyE_IgG", "HlyE_IgA"), + lambda = 0.1 + ) %>% + print() -print(ll_A + ll_G) +print(ll_a + ll_g) ``` ## graph log-likelihood -We can also graph the log-likelihoods, even without finding the MLEs, using `graph.loglik()`: +We can also graph the log-likelihoods, even without finding the MLEs, using `graph_loglik()`: ```{r} -lik_HlyE_IgA <- graph.loglik( - pop_data = csdata, - curve_params = dmcmc, - noise_params = cond, - antigen_isos = "HlyE_IgA", - log_x = TRUE -) +lik_HlyE_IgA <- + graph_loglik( + pop_data = csdata, + curve_params = dmcmc, + noise_params = cond, + antigen_isos = "HlyE_IgA", + log_x = TRUE + ) -lik_HlyE_IgG <- graph.loglik( +lik_HlyE_IgG <- graph_loglik( previous_plot = lik_HlyE_IgA, pop_data = csdata, curve_params = dmcmc, @@ -216,7 +237,7 @@ lik_HlyE_IgG <- graph.loglik( log_x = TRUE ) -lik_both <- graph.loglik( +lik_both <- graph_loglik( previous_plot = lik_HlyE_IgG, pop_data = csdata, curve_params = dmcmc, @@ -238,9 +259,10 @@ est1 <- est.incidence( curve_params = dmcmc, noise_params = cond, lambda_start = .1, - build_graph = T, - verbose = T, # print updates as the function runs - print_graph = F, # display the log-likelihood curve while `est.incidence()` is running + build_graph = TRUE, + verbose = TRUE, # print updates as the function runs + print_graph = FALSE, # display the log-likelihood curve while + #`est.incidence()` is running antigen_isos = antibodies ) ``` @@ -268,7 +290,6 @@ autoplot(est1, log_x = TRUE) ```{r "init-parallel"} library(parallel) n_cores <- max(1, parallel::detectCores() - 1) -# n_cores = 1 print(n_cores) ``` @@ -281,13 +302,12 @@ nclus <- 10 nrep <- 100 # incidence rate in e -lmbdaVec <- c(.05, .1, .15, .2) +lambdas <- c(.05, .1, .15, .2) -sim.df <- +sim_df <- sim.cs.multi( - # verbose = TRUE, n_cores = n_cores, - lambdas = lmbdaVec, + lambdas = lambdas, nclus = nclus, n.smpl = nrep, age.rng = lifespan, @@ -299,16 +319,17 @@ sim.df <- format = "long" ) -print(sim.df) +print(sim_df) ``` We can plot the distributions of the simulated responses: ```{r} -ggplot(sim.df, aes( - x = as.factor(cluster), - y = value -)) + +sim_df %>% ggplot() + + aes( + x = as.factor(cluster), + y = value + ) + geom_beeswarm(size = .2, alpha = .3, aes(color = antigen_iso)) + geom_boxplot(outlier.colour = NA, fill = NA) + scale_y_log10() + @@ -320,9 +341,8 @@ ggplot(sim.df, aes( ```{r, "est-by-stratum"} ests <- - # sim.df %>% est.incidence.by( - pop_data = sim.df, + pop_data = sim_df, curve_params = dmcmc, noise_params = cond, num_cores = n_cores, @@ -332,7 +352,6 @@ ests <- verbose = TRUE, build_graph = TRUE, # slows down the function substantially antigen_isos = c("HlyE_IgG", "HlyE_IgA") - # antigen_isos = "HlyE_IgA" ) ``` diff --git a/vignettes/articles/templates/serocalculator_GatesIDM_10.01.24.potx b/vignettes/articles/templates/serocalculator_GatesIDM_10.01.24.potx new file mode 100644 index 00000000..5a797572 Binary files /dev/null and b/vignettes/articles/templates/serocalculator_GatesIDM_10.01.24.potx differ diff --git a/vignettes/fig/fig1a-1.svg b/vignettes/fig/fig1a-1.svg new file mode 100644 index 00000000..f68e997e --- /dev/null +++ b/vignettes/fig/fig1a-1.svg @@ -0,0 +1,3795 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/vignettes/fig/response.pdf b/vignettes/fig/response.pdf new file mode 100644 index 00000000..f4a78eb2 Binary files /dev/null and b/vignettes/fig/response.pdf differ diff --git a/vignettes/man/figures/qr_akm2.svg b/vignettes/man/figures/qr_akm2.svg new file mode 100644 index 00000000..2e6e89cc --- /dev/null +++ b/vignettes/man/figures/qr_akm2.svg @@ -0,0 +1,535 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/vignettes/methodology.qmd b/vignettes/methodology.qmd new file mode 100644 index 00000000..b9cda730 --- /dev/null +++ b/vignettes/methodology.qmd @@ -0,0 +1,353 @@ +--- +title: "Introduction to seroincidence estimation" +vignette: > + %\VignetteIndexEntry{Introduction to seroincidence estimation} + %\VignetteEngine{quarto::html} + %\VignetteEncoding{UTF-8} +--- + +{{< include articles/_macros.qmd >}} + +# Overview + +## Defining incidence + +The **incidence rate** of a disease over a **specific time period** is +the rate at which individuals in a population are *acquiring* the disease +during that time period [@Noordzij2010diseasemeasures]. + +**Example:** +if there are **10 new cases** of typhoid +in a population of **1000 persons** +during a **one month** time period, +then the incidence rate for that time period is +**10 new cases per 1000 persons per month**. + +## Mathematical definition of incidence {visibility="hidden"} + +More mathematically, +the incidence rate *at a given time point* is +the *derivative* (i.e., the current rate of change over time) +of the expected cumulative count of infections per person at risk, +at that time: + +$$\frac{d}{dt} \mathbb{E}\left[\frac{C(t)}{n} \mid N(t) =n\right]$$ + +where $C(t)$ is the cumulative total number of infections in +the population of interest, +and $N(t)$ is the number of individuals at risk at time $t$. + +## Scale of incidence rates {visibility="hidden"} + +In both definitions, +the units for an incidence rate are +"# new infections per # persons at risk per time duration"; +for example, "new infections per 1000 persons per year". + +For convenience, +we can rescale the incidence rate to make it easier to understand; +for example, we might express incidence as +"# new infections per 1000 persons per year" +or "# new infections per 100,000 persons per day", etc. + +## Incidence from an individual's perspective + +From the perspective of an individual in the population: + +* the **incidence rate** (at a given time point ($t$) +is the instantaneous **probability** (density) of +**becoming infected** at that time point, +**given** that they are **at risk** at that time point. + +* That is, the incidence rate is a **hazard** rate. + +* Notation: let's use **$\lambda_{t}$** to denote the incidence rate at time $t$. + +# Estimating incidence from cross-sectional antibody surveys + +## Cross-sectional antibody surveys {.incremental} + +::: notes +Typically, it is difficult to estimate changes from a single time point. +However, we can sometimes make assumptions that allow us to do so. +In particular, if we assume that the incidence rate is constant over time, +then we can estimate incidence from a single cross-sectional survey. + +We will need two pieces of notation to formalize this process. + +::: + +::: incremental + +* We recruit participants from the population of interest. + +* For each survey participant, we measure antibody levels $(Y)$ +for the disease of interest + +* Each participant was **most recently infected** at some time $(T)$ +**prior** to when we measured their antibodies. + +::: notes +* If a participant has never been infected since birth, then $T$ is undefined. + + +::: + +* $T$ is a **latent, unobserved variable**. + +::: notes +* We **don't directly observe $T$**; we **only observe $Y$**, +which we hope tells us something about $T$ and $\lambda$. + +::: + +::: + +## Modeling assumptions + +We **assume** that: + +::: incremental + +* The incidence rate is approximately **constant** +**over time** and **across the population** ("**constant and homogenous incidence**") + +::: notes + +* that is: +$$\lambda_{i,t} = \lambda, \forall i,t$$ + +(We can analyze subpopulations separately to make homogeneity more plausible.) +::: + +* Participants are always at risk of a new infection, +regardless of how recently they have been infected +("**no lasting immunity**"). + +::: + +::: notes +(For diseases like typhoid, +the no-immunity assumption may not hold exactly, +but hopefully approximately; +modeling the effects of re-exposure during an active infection is +[on our to-do list](https://github.com/UCD-SERG/dcm/issues/11)). +::: + +## Time since infection and incidence + +Under those assumptions: + +::: incremental + +* $T$ has an **exponential distribution**: + +$$\pdf(T=t) = \lambda \expf{-\lambda t}$$ + + * More precisely, the distribution is exponential + **truncated by age** at observation ($a$): + +$$\pdf(T=t|A=a) = 1_{t \in[0,a]}\lambda \expf{-\lambda t} + 1_{t = \NA} \expf{-\lambda a}$$ + +* the rate parameter $\lambda$ is the incidence rate + +::: + +::: notes +This is a time-to-event model, +looking **backwards in time** from the survey date +(when the blood sample was collected). + +The probability that +an individual was **last** infected $t$ days ago, $p(T=t)$, +is equal to the probability of being infected at time $t$ +(i.e., the incidence rate at time $t$, $\lambda$) +times the probability of not being infected after time $t$, +which turns out to be $\exp(-\lambda t)$. + +The distribution of $T$ is truncated by the patient's birth date; +the probability that they have never been infected is +$\expf{-\lambda a}$, +where $a$ is the patient's age at the time of the survey. +::: + +## Likelihood of latent infection times {.smaller} + +{{< include articles/_sec-latent-likelihood.qmd >}} + +## Example log-likelihood curves + +::: notes +Here's what that would look like: +::: + +```{r} +#| fig-cap: "Example log-likelihood curves" +#| label: fig-ex-lik-curves + +library(serocalculator) +library(dplyr) +# Import longitudinal antibody parameters from OSF +curves <- + "https://osf.io/download/rtw5k/" %>% + load_curve_params() %>% + filter(iter < 50) + +# Import cross-sectional data from OSF and rename required variables: +xs_data <- + "https://osf.io/download//n6cp3/" %>% + load_pop_data() + +noise <- url("https://osf.io/download//hqy4v/") %>% readRDS() +``` + +```{r} +#| fig-cap: "Example log(likelihood) curves" +#| label: fig-loglik +lik_HlyE_IgA <- graph_loglik( + pop_data = xs_data, + curve_params = curves, + noise_params = noise, + antigen_isos = "HlyE_IgA", + log_x = TRUE +) + +lik_HlyE_IgG <- graph_loglik( + previous_plot = lik_HlyE_IgA, + pop_data = xs_data, + curve_params = curves, + noise_params = noise, + antigen_isos = "HlyE_IgG", + log_x = TRUE +) + +lik_both <- graph_loglik( + previous_plot = lik_HlyE_IgG, + pop_data = xs_data, + curve_params = curves, + noise_params = noise, + antigen_isos = c("HlyE_IgG", "HlyE_IgA"), + log_x = TRUE +) + +print(lik_both) +``` + +## standard error + +The standard error of the estimate is approximately equal to the inverse of the rate of curvature (2nd derivative, aka Hessian) in the log-likelihood function, at the maximum: + +more curvature -> likelihood peak is clearer -> smaller standard errors + +## Likelihood of observed data + +::: notes +Unfortunately, we don't observe infection times $T$; +we only observe antibody levels ${Y}$. +So things get a little more complicated. + +In short, we are hoping that we can estimate +$T$ (time since last infection) +from $Y$ (current antibody levels). +If we could do that, +then we could plug in our estimates $\hat t_i$ +into that likelihood above, +and estimate $\lambda$ as previously. + +We're actually going to do something a little more nuanced; +instead of just using one value for $\hat t$, +we are going to consider all possible values of $t$ for each individual. + +We need to link the data we actually observed to the incidence rate. + +The likelihood of an individual's observed data, +$\pdf(Y=y)$, +can be expressed as an integral over +the joint likelihood of $Y$ and $T$ +(using the Law of Total Probability): +::: + +::: incremental + +* $$\pdf(Y=y) = \int_t \pdf(Y=y,T=t)dt$$ + +:::: notes +Further, we can express the joint probability $p(Y=y,T=t)$ as the product of +$p(T=t)$ and +$p(Y=y|T=t)$ the "antibody response curve after infection". +That is: +:::: + +* $$\pdf(Y=y,T=t) = \pdf(Y=y|T=t) \pdf(T=t)$$ + +::: + +## Antibody response curves + +::: {#fig-decay} + +![](fig/fig1a-1.svg) + +Antibody response curves, $p(Y=y|T=t)$, for typhoid + +::: + +## Putting it all together + +Substituting $p(Y=y,T=t) = p(Y=y|T=t)P(T=t)$ +into the previous expression for $p(Y=y)$: + +$$ +\begin{aligned} +p(Y=y) +&= \int_t p(Y=y|T=t)P(T=t) dt +\end{aligned} +$$ + +## Composing the likelihood {.smaller} + +::: notes +Now, the likelihood of the observed data +$\vec{y} = (y_1, y_2, ..., y_n)$ is: +::: + +$$ +\begin{aligned} +\mathcal{L}(\lambda) +&= \prod_{i=1}^n p(Y=y_i) +\\&= \prod_{i=1}^n \int_t p(Y=y_i|T=t)p_\lambda(T=t)dt\\ +\end{aligned} +$$ + +::: notes +If we know $p(Y=y|T=t)$, then we can maximize $\mathcal{L}(\lambda)$ over $\lambda$ to find the "maximum likelihood estimate" (MLE) of $\lambda$, denoted $\hat\lambda$. +::: + +## Finding the MLE numerically + +::: notes +The likelihood of $Y$ involves the product of integrals, so the log-likelihood involves the sum of the logs of integrals: +::: + +$$ +\begin{aligned} +\log \mathcal{L} (\lambda) +&= \log \prod_{i=1}^n \int_t p(Y=y_i|T=t)p_\lambda(T=t)dt\\ +&= \sum_{i=1}^n \log\left\{\int_t p(Y=y_i|T=t)p_\lambda(T=t)dt\right\}\\ +\end{aligned} +$$ + +::: notes + +The derivative of this expression doesn't come out cleanly, so we will use a *numerical method* (specifically, a Newton-type algorithm, implemented by `stats::nlm()`) to find the MLE and corresponding standard error. + +::: + +{{< include articles/_antibody-response-model.qmd >}} + + + +## References {.smaller .unnumbered} + +::: {#refs} +::: diff --git a/vignettes/references.bib b/vignettes/references.bib new file mode 100644 index 00000000..7320f3d6 --- /dev/null +++ b/vignettes/references.bib @@ -0,0 +1,173 @@ +@article{de_Graaf_2014, + doi = {10.1016/j.epidem.2014.08.002}, + url = {https://doi.org/10.1016%2Fj.epidem.2014.08.002}, + year = 2014, + month = {dec}, + publisher = {Elsevier {BV}}, + volume = {9}, + pages = {1--7}, + author = {W.F. deGraaf and M.E.E. Kretzschmar and P.F.M. Teunis and O. Diekmann}, + title = {A two-phase within-host model for immune response and its application to serological profiles of pertussis}, + journal = {Epidemics} +} + +@article{Simonsen_2009, + doi = {10.1002/sim.3592}, + url = {https://doi.org/10.1002%2Fsim.3592}, + year = 2009, + month = {jun}, + publisher = {Wiley}, + volume = {28}, + number = {14}, + pages = {1882--1895}, + author = {J. Simonsen and K. M{\o}lbak and G. Falkenhorst and K. A. Krogfelt and A. Linneberg and P.F.M. Teunis}, + title = {Estimation of incidences of infectious diseases based on antibody measurements}, + journal = {Statistics in Medicine} +} + +@article{Teunis_2012, + doi = {10.1002/sim.5322}, + url = {https://doi.org/10.1002%2Fsim.5322}, + year = 2012, + month = {mar}, + publisher = {Wiley}, + volume = {31}, + number = {20}, + pages = {2240--2248}, + author = {P.F.M. Teunis and JCH van Eijkeren and CW Ang and YTHP van Duynhoven and JB Simonsen and MA Strid and W van Pelt}, + title = {Biomarker dynamics: estimating infection rates from serological data}, + journal = {Statistics in Medicine} +} + +@article{Teunis_2016, + doi = {10.1016/j.epidem.2016.04.001}, + url = {https://doi.org/10.1016%2Fj.epidem.2016.04.001}, + year = 2016, + month = {sep}, + publisher = {Elsevier {BV}}, + volume = {16}, + pages = {33--39}, + author = {P.F.M. Teunis and J.C.H. van Eijkeren and W.F. de Graaf and A. Bona{\v{c}}i{\'{c}} Marinovi{\'{c}} and M.E.E. Kretzschmar}, + title = {Linking the seroresponse to infection to within-host heterogeneity in antibody production}, + journal = {Epidemics} +} + +@article{Strid_2001, + doi = {10.1128/cdli.8.2.314-319.2001}, + url = {https://doi.org/10.1128%2Fcdli.8.2.314-319.2001}, + year = 2001, + month = {mar}, + publisher = {American Society for Microbiology}, + volume = {8}, + number = {2}, + pages = {314--319}, + author = {Mette Aagaard Strid and J{\o}rgen Engberg and Lena Brandt Larsen and Kamilla Begtrup and K{\aa}re M{\o}lbak and Karen Angeliki Krogfelt}, + title = {Antibody responses to Campylobacter infections determined by an enzyme-linked immunosorbent assay: 2-year follow-up study of 210 patients}, + journal = {Clinical Diagnostic Laboratory Immunology} +} + +@article{satten_etal04, + author = {Satten, Glen A. and Kong, Fanhui and Wright, David J. and Glynn, Simone A. and Schreiber, George B.}, + title = "{How special is a ‘special’ interval: modeling departure from length‐biased sampling in renewal processes}", + journal = {Biostatistics}, + volume = {5}, + number = {1}, + pages = {145-151}, + year = {2004}, + month = {01}, + issn = {1465-4644}, + doi = {10.1093/biostatistics/5.1.145}, + url = {https://doi.org/10.1093/biostatistics/5.1.145}, + eprint = {https://academic.oup.com/biostatistics/article-pdf/5/1/145/941872/050145.pdf}, +} + +@article{zelen04, + title={Forward and backward recurrence times and length biased sampling: age specific models}, + author={Zelen, Marvin}, + journal={Lifetime Data Analysis}, + volume={10}, + pages={325--334}, + year={2004}, + publisher={Springer} +} + + +@article{degreeff_etal12, + title={Two-component cluster analysis of a large serodiagnostic database for specificity of increases of IgG antibodies against pertussis toxin in paired serum samples and of absolute values in single serum samples}, + author={de Greeff, Sabine C and Teunis, Peter and de Melker, Hester E and Mooi, Frits R and Notermans, Daan W and Elvers, Bert and Schellekens, Joop FP}, + journal={Clinical and Vaccine Immunology}, + volume={19}, + number={9}, + pages={1452--1456}, + year={2012}, + publisher={Am Soc Microbiol}} + +@article{Aiemjoy_2022_Lancet, + url = {https://doi.org/10.1016/S2666-5247(22)00114-8}, + year = 2022, + volume = {3}, + number = {8}, + pages = {e578--e587}, + author = {K. Aiemjoy and Seidman J. C. and Saha S. and Munira S. J. and Islam Sajib M. S. and Sium, S. M. al, Sarkar, A., Alam, N., Zahan, F. N., Kabir, M. S., Tamrakar, D., Vaidya, K., Shrestha, R., Shakya, J., Katuwal, N., Shrestha, S., Yousafzai, M. T., Iqbal, J., Dehraj, I. F., … Andrews, J. R.}, + title = {Estimating typhoid incidence from community-based serosurveys: a multicohort study}, + journal = {The Lancet Microbe} +} + +@article{Aiemjoy_2024_scrub, + title={Estimating the seroincidence of scrub typhus using antibody dynamics following infection}, + author={Aiemjoy, Kristen and Katuwal, Nishan and Vaidya, Krista and Shrestha, Sony and Thapa, Melina and Teunis, Peter and Bogoch, Isaac I and Trowbridge, Paul and Kantipong, Pacharee and Blacksell, Stuart D and others}, + journal={American Journal of Tropical Medicine and Hygiene}, + year={Accepted Feb 2024}, + +} + +@article{Aiemjoy_2022_SouthSudan, + doi = {10.3201/eid2811.220239}, + year = 2022, + month = {nov}, + volume = {28}, + number = {11}, + pages = {2316--2320}, + author = {K. Aiemjoy and John Rumunu and Juma John Hassen, Kirsten E. Wiens, Denise Garrett, Polina Kamenskaya, Jason B. Harris, Andrew S. Azman, Peter Teunis, +Jessica C. Seidman, Joseph F. Wamala, Jason R. Andrews, Richelle C. Charles}, + title = {Seroincidence of Enteric Fever,Juba, South Sudan}, + journal = {Emerging Infectious Diseases} +} + +@article{Teunis_2020, + doi = {10.1002/sim.8578}, + year = 2020, + volume = {39}, + number = {21}, + pages = {2799--2814}, + author = { P.F.M. Teunis and J.C.H. van Eijkeren}, + title = {Estimation of seroconversion rates for infectious diseases: Effects of age and noise}, + journal = {Statistics in Medicine} +} + +@article{Noordzij2010diseasemeasures, + author = {Noordzij, Marlies and Dekker, Friedo W. and Zoccali, Carmine and Jager, Kitty J.}, + title = "{Measures of Disease Frequency: Prevalence and Incidence}", + journal = {Nephron Clinical Practice}, + volume = {115}, + number = {1}, + pages = {c17-c20}, + year = {2010}, + month = {02}, + abstract = "{To describe how often a disease or another health event occurs in a population, different measures of disease frequency can be used. The prevalence reflects the number of existing cases of a disease. In contrast to the prevalence, the incidence reflects the number of new cases of disease and can be reported as a risk or as an incidence rate. Prevalence and incidence are used for different purposes and to answer different research questions. In this article, we discuss the different measures of disease frequency and we explain when to apply which measure.}", + issn = {1660-2110}, + doi = {10.1159/000286345}, + url = {https://doi.org/10.1159/000286345}, + eprint = {https://karger.com/nec/article-pdf/115/1/c17/3766014/000286345.pdf}, +} + +@article{volterra1928variations, + title={Variations and fluctuations of the number of individuals in animal species living together}, + author={Volterra, Vito}, + journal={ICES Journal of Marine Science}, + volume={3}, + number={1}, + pages={3--51}, + year={1928}, + publisher={Oxford University Press} +}