serocalculator: Estimating Infection Rates from Serological Data

Description

Translates antibody levels measured in cross-sectional population samples into estimates of the frequency with which seroconversions (infections) occur in the sampled populations. Replaces the previous 'seroincidence' package.

Author(s)

Maintainer: Kristina Lai kwlai@ucdavis.edu

Authors:

• Peter Teunis p.teunis@emory.edu (Author of the method and original code.) [copyright holder]

• Chris Orwa

• Kristen Aiemjoy kaiemjoy@ucdavis.edu

• Douglas Ezra Morrison demorrison@ucdavis.edu

Other contributors:

• Kwan Ho Lee ksjlee@ucdavis.edu [contributor]

See Also

Useful links:

• https://ucd-serg.github.io/serocalculator/

• https://github.com/UCD-SERG/serocalculator

• Report bugs at https://github.com/UCD-SERG/serocalculator/issues

Calculate negative log-likelihood

Description

Same as log_likelihood(), except negated and requiring lambda on log scale (used in combination with nlm(), to ensure that the optimization search doesn't stray into negative values of lambda).

Usage

.nll(log.lambda, ...)

Arguments

Value

the negative log-likelihood of the data with the current parameter values

Extract or replace parts of a seroincidence.by object

Description

Extract or replace parts of a seroincidence.by object

Usage

## S3 method for class 'seroincidence.by'
x[i, ...]

Arguments

Value

the subset specified

kinetics of the antibody (ab) response (power function decay)

Description

kinetics of the antibody (ab) response (power function decay)

Usage

ab(t, par, ...)

Arguments

Value

a matrix() of predicted biomarker values

Examples

par1 <- matrix(
    c(
      1.11418923843475, 1, 0.12415057798022207, 0.24829344792968783,
      0.01998946878312856, 0.0012360802436587237, 1.297194045996013,
      1.3976510415108334, 1, 0.2159993563893431, 0.4318070551383313,
      0.0015146395107173347, 0.0003580062906750277, 1.5695811573082081
    ),
    nrow = 7L,
    ncol = 2L,
    dimnames = list(
      params = c("y0", "b0", "mu0", "mu1", "c1", "alpha", "shape_r"),
      antigen_iso = c("HlyE_IgA", "HlyE_IgG")
    )
    )
t <- 0:1444
blims <- matrix(
   rep(c(0, 0.5), each = 2L),
   nrow = 2L,
   ncol = 2L,
   dimnames = list(c("HlyE_IgA", "HlyE_IgG"), c("min", "max"))
   )
preds <- serocalculator:::ab(t = t, par = par1, blims = blims)

Analyze simulation results

Description

Analyze simulation results

Usage

analyze_sims(data)

Arguments

Value

a sim_results object (extends tibble::tbl_df)

Examples

dmcmc <- typhoid_curves_nostrat_100

n_cores <- 2

nclus <- 20
# cross-sectional sample size
nrep <- c(50, 200)

# incidence rate in e
lambdas <- c(.05, .8)

antibodies <- c("HlyE_IgA", "HlyE_IgG")
lifespan <- c(0, 10)
dlims = rbind(
"HlyE_IgA" = c(min = 0, max = 0.5),
"HlyE_IgG" = c(min = 0, max = 0.5)
)
sim_df <- sim_pop_data_multi(
n_cores = n_cores,
lambdas = lambdas,
nclus = nclus,
sample_sizes = nrep,
age_range = lifespan,
antigen_isos = antibodies,
renew_params = FALSE,
add_noise = TRUE,
curve_params = dmcmc,
noise_limits = dlims,
format = "long"
)
cond <- tibble::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)
)
ests <-
est_seroincidence_by(
pop_data = sim_df,
sr_params = dmcmc,
noise_params = cond,
num_cores = n_cores,
strata = c("lambda.sim", "sample_size", "cluster"),
curve_strata_varnames = NULL,
noise_strata_varnames = NULL,
verbose = FALSE,
build_graph = FALSE, # slows down the function substantially
antigen_isos = c("HlyE_IgG", "HlyE_IgA")
)

ests |>
summary() |>
analyze_sims()

Load antibody decay curve parameter

Description

as_curve_params() was renamed to as_sr_params() to create a more consistent API.

Usage

as_curve_params(...)

Load noise parameters

Description

Load noise parameters

Usage

as_noise_params(data, antigen_isos = NULL)

Arguments

Value

a noise_params object (a tibble::tbl_df with extra attribute antigen_isos)

Examples

library(magrittr)
noise_data <-
  serocalculator_example("example_noise_params.csv") %>%
  read.csv() %>%
  as_noise_params()

print(noise_data)

Load a cross-sectional antibody survey data set

Description

Load a cross-sectional antibody survey data set

Usage

as_pop_data(
  data,
  antigen_isos = NULL,
  age = "Age",
  value = "result",
  id = "index_id",
  standardize = TRUE
)

Arguments

Value

a pop_data object (a tibble::tbl_df with extra attribute antigen_isos)

Examples

library(magrittr)
xs_data <-
  serocalculator_example("example_pop_data.csv") |>
  read.csv() |>
  as_pop_data()

print(xs_data)

Load longitudinal seroresponse parameters

Description

Load longitudinal seroresponse parameters

Usage

as_sr_params(data, antigen_isos = NULL)

Arguments

Value

a curve_data object (a tibble::tbl_df with extra attribute antigen_isos)

Examples

library(magrittr)
curve_data <-
  serocalculator_example("example_curve_params.csv") %>%
  read.csv() %>%
  as_curve_params()

print(curve_data)

Graph antibody decay curves by antigen isotype

Description

Graph antibody decay curves by antigen isotype

Usage

## S3 method for class 'curve_params'
autoplot(
  object,
  method = c("graph.curve.params", "graph_seroresponse_model_1"),
  ...
)

Arguments

Details

Currently, the backend for this method is graph.curve.params(). Previously, the backend for this method was graph_seroresponse_model_1(). That function is still available if preferred.

Value

a ggplot2::ggplot() object

Examples

library(dplyr)
library(ggplot2)
library(magrittr)

curve <-
  serocalculator_example("example_curve_params.csv") |>
  read.csv() |>
  as_sr_params() |>
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG")) |>
  autoplot()

curve

Plot distribution of antibodies

Description

autoplot() method for pop_data objects

Usage

## S3 method for class 'pop_data'
autoplot(object, log = FALSE, type = "density", strata = NULL, ...)

Arguments

Value

a ggplot2::ggplot object

Examples

library(dplyr)
library(ggplot2)
library(magrittr)

xs_data <-
  serocalculator_example("example_pop_data.csv") |>
  read.csv() |>
  as_pop_data()

xs_data |> autoplot(strata = "catchment", type = "density")
xs_data |> autoplot(strata = "catchment", type = "age-scatter")

Plot the log-likelihood curve for the incidence rate estimate

Description

Plot the log-likelihood curve for the incidence rate estimate

Usage

## S3 method for class 'seroincidence'
autoplot(object, log_x = FALSE, ...)

Arguments

Value

a ggplot2::ggplot()

Examples

library(dplyr)
library(ggplot2)

xs_data <-
  sees_pop_data_pk_100

curve <-
  typhoid_curves_nostrat_100 %>%
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

noise <-
  example_noise_params_pk

est1 <- est_seroincidence(
  pop_data = xs_data,
  sr_param = curve,
  noise_param = noise,
  antigen_isos = c("HlyE_IgG", "HlyE_IgA"),
  build_graph = TRUE
)

# Plot the log-likelihood curve
autoplot(est1)

Plot seroincidence.by log-likelihoods

Description

Plots log-likelihood curves by stratum, for seroincidence.by objects

Usage

## S3 method for class 'seroincidence.by'
autoplot(object, ncol = min(3, length(object)), ...)

Arguments

Value

a "ggarrange" object: a single or list() of ggplot2::ggplot()s

Examples

library(dplyr)
library(ggplot2)

xs_data <-
  sees_pop_data_pk_100

curve <-
  typhoid_curves_nostrat_100 %>%
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

noise <-
  example_noise_params_pk

est2 <- est_seroincidence_by(
  strata = c("catchment"),
  pop_data = xs_data,
  sr_params = curve,
  curve_strata_varnames= NULL,
  noise_strata_varnames = NULL,
  noise_params = noise,
  antigen_isos = c("HlyE_IgG", "HlyE_IgA"),
  #num_cores = 8, #Allow for parallel processing to decrease run time
  build_graph = TRUE
)

# Plot the log-likelihood curve
autoplot(est2)

Plot simulation results autoplot() method for sim_results objects

Description

Plot simulation results autoplot() method for sim_results objects

Usage

## S3 method for class 'sim_results'
autoplot(object, statistic = "Empirical_SE", ...)

Arguments

Value

a ggplot2::ggplot

Examples

dmcmc <- typhoid_curves_nostrat_100

n_cores <- 2

nclus <- 20
# cross-sectional sample size
nrep <- c(50, 200)

# incidence rate in e
lambdas <- c(.05, .8)
lifespan <- c(0, 10)
antibodies <- c("HlyE_IgA", "HlyE_IgG")
dlims <- rbind(
"HlyE_IgA" = c(min = 0, max = 0.5),
"HlyE_IgG" = c(min = 0, max = 0.5)
)
sim_df <-
sim_pop_data_multi(
n_cores = n_cores,
lambdas = lambdas,
nclus = nclus,
sample_sizes = nrep,
age_range = lifespan,
antigen_isos = antibodies,
renew_params = FALSE,
add_noise = TRUE,
curve_params = dmcmc,
noise_limits = dlims,
format = "long"
)
cond <- tibble::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)
)
ests <-
est_seroincidence_by(
pop_data = sim_df,
sr_params = dmcmc,
noise_params = cond,
num_cores = n_cores,
strata = c("lambda.sim", "sample_size", "cluster"),
curve_strata_varnames = NULL,
noise_strata_varnames = NULL,
verbose = FALSE,
build_graph = FALSE, # slows down the function substantially
antigen_isos = c("HlyE_IgG", "HlyE_IgA")
)

ests |>
summary() |>
analyze_sims() |>
autoplot()

Plot method for summary.seroincidence.by objects

Description

Plot method for summary.seroincidence.by objects

Usage

## S3 method for class 'summary.seroincidence.by'
autoplot(object, type, ...)

Arguments

Value

a ggplot2::ggplot() object

Examples

library(dplyr)
library(ggplot2)

xs_data <-
  sees_pop_data_pk_100

curve <-
  typhoid_curves_nostrat_100 %>%
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

noise <-
  example_noise_params_pk

est2 <- est_seroincidence_by(
  strata = c("catchment", "ageCat"),
  pop_data = xs_data,
  sr_params = curve,
  noise_params = noise,
  curve_strata_varnames= NULL,
  noise_strata_varnames = NULL,
  antigen_isos = c("HlyE_IgG", "HlyE_IgA"),
  num_cores = 2 # Allow for parallel processing to decrease run time
)

est2sum <- summary(est2)

est2sum |> autoplot(
    type ="scatter",
    xvar = "ageCat",
    color_var = "catchment",
    CIs = TRUE,
    group_var = "catchment")

est2sum |> autoplot(
    type = "bar",
    yvar = "ageCat",
    color_var = "catchment",
    CIs = TRUE)

Substitute baseline values

Description

whenever y is below a cutoff (blims[kab,2]), substitute a random sample from a baseline distribution

Usage

baseline(kab, yvec, blims, ...)

Arguments

Value

an altered version of yvec

Check the formatting of a cross-sectional antibody survey dataset.

Description

Check the formatting of a cross-sectional antibody survey dataset.

Usage

check_pop_data(pop_data, verbose = FALSE)

Arguments

Value

NULL (invisibly)

Examples

library(magrittr)

xs_data <-
  serocalculator_example("example_pop_data.csv") |>
  read.csv() |>
  as_pop_data()

check_pop_data(xs_data, verbose = TRUE)

Count observations by stratum

Description

Count observations by stratum

Usage

count_strata(
  data,
  strata_varnames,
  biomarker_names_var = get_biomarker_names_var(data)
)

Arguments

Value

a tibble::tbl_df counting observations by stratum

Examples

sees_pop_data_pk_100 |> count_strata(strata_varnames = "catchment")

Convert a data.frame (or tibble) into a multidimensional array

Description

df.to.array() was renamed to df_to_array() to create a more consistent API.

Usage

df.to.array(df, dim_var_names, value_var_name = "value")

Convert a data.frame (or tibble) into a multidimensional array

Description

Convert a data.frame (or tibble) into a multidimensional array

Usage

df_to_array(df, dim_var_names, value_var_name = "value")

Arguments

Value

an array() with dimensions defined by the variables in df listed in dim_var_names

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])

Estimate Seroincidence

Description

est.incidence() was renamed to est_seroincidence() to create a more consistent API.

Usage

est.incidence(...)

Estimate Seroincidence

Description

est.incidence.by() was renamed to est_seroincidence_by() to create a more consistent API.

Usage

est.incidence.by(...)

Find the maximum likelihood estimate of the incidence rate parameter

Description

This function models seroincidence using maximum likelihood estimation; that is, it finds the value of the seroincidence parameter which maximizes the likelihood (i.e., joint probability) of the data.

Usage

est_seroincidence(
  pop_data,
  sr_params,
  noise_params,
  antigen_isos = get_biomarker_names(pop_data),
  lambda_start = 0.1,
  stepmin = 1e-08,
  stepmax = 3,
  verbose = FALSE,
  build_graph = FALSE,
  print_graph = build_graph & verbose,
  ...
)

Arguments

Value

a "seroincidence" object, which is a stats::nlm() fit object with extra metadata attributes lambda_start, antigen_isos, and ll_graph

Examples

library(dplyr)

xs_data <-
  sees_pop_data_pk_100

sr_curve <-
  typhoid_curves_nostrat_100 |>
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

noise <-
  example_noise_params_pk

est1 <- est_seroincidence(
  pop_data = xs_data,
  sr_params = sr_curve,
  noise_params = noise,
  antigen_isos = c("HlyE_IgG", "HlyE_IgA"),
)

summary(est1)

Estimate Seroincidence

Description

Function to estimate seroincidences based on cross-sectional serology data and longitudinal response model.

Usage

est_seroincidence_by(
  pop_data,
  sr_params,
  noise_params,
  strata,
  curve_strata_varnames = strata,
  noise_strata_varnames = strata,
  antigen_isos = unique(pull(pop_data, "antigen_iso")),
  lambda_start = 0.1,
  build_graph = FALSE,
  num_cores = 1L,
  verbose = FALSE,
  print_graph = FALSE,
  ...
)

Arguments

Details

If strata is left empty, a warning will be produced, recommending that you use est_seroincidence() for unstratified analyses, and then the data will be passed to est_seroincidence(). If for some reason you want to use est_seroincidence_by() with no strata instead of calling est_seroincidence(), you may use NA, NULL, or "" as the strata argument to avoid that warning.

Value

• if strata has meaningful inputs: An object of class "seroincidence.by"; i.e., a list of "seroincidence" objects from est_seroincidence(), one for each stratum, with some meta-data attributes.

• if strata is missing, NULL, NA, or "": An object of class "seroincidence".

Examples

library(dplyr)

xs_data <-
  sees_pop_data_pk_100

curve <-
  typhoid_curves_nostrat_100 |>
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

noise <-
  example_noise_params_pk

est2 <- est_seroincidence_by(
  strata = "catchment",
  pop_data = xs_data,
  sr_params = curve,
  noise_params = noise,
  antigen_isos = c("HlyE_IgG", "HlyE_IgA"),
  # num_cores = 8 # Allow for parallel processing to decrease run time
  iterlim = 5 # limit iterations for the purpose of this example
)
print(est2)
summary(est2)

Small example of noise parameters for typhoid

Description

A subset of noise parameter estimates from the SEES study, for examples and testing, for Pakistan

Usage

example_noise_params_pk

Format

example_noise_params_pk

A curve_params object (from as_sr_params()) with 4 rows and 7 columns:

antigen_iso

which antigen and isotype are being measured (data is in long format)

Country

Location for which the noise parameters were estimated

y.low

Lower limit of detection

eps

Measurement noise, 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.

nu

Biological noise: error from cross-reactivity to other antibodies. It is defined as the 95th percentile of the distribution of antibody responses to the antigen-isotype in a population with no exposure.

y.high

Upper limit of detection

Lab

Lab for which noise was estimated.

Source

https://osf.io/rtw5k

Small example of noise parameters for typhoid

Description

A subset of noise parameter estimates from the SEES study, for examples and testing.

Usage

example_noise_params_sees

Format

example_noise_params_pk

A curve_params object (from as_sr_params()) with 4 rows and 7 columns:

antigen_iso

which antigen and isotype are being measured (data is in long format)

Country

Location for which the noise parameters were estimated

y.low

Lower limit of detection

eps

Measurement noise, 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.

nu

Biological noise: error from cross-reactivity to other antibodies. It is defined as the 95th percentile of the distribution of antibody responses to the antigen-isotype in a population with no exposure.

y.high

Upper limit of detection

Lab

Lab for which noise was estimated.

Source

https://osf.io/rtw5k

Snapshot testing for data.frames

Description

copied from https://github.com/bcgov/ssdtools with permission (https://github.com/bcgov/ssdtools/issues/379)

Usage

expect_snapshot_data(x, name, digits = 6)

Arguments

Value

NULL (from testthat::expect_snapshot_file())

Examples

## Not run: 
expect_snapshot_data(iris, name = "iris")

## End(Not run)

Calculate negative log-likelihood (deviance) for one antigen:isotype pair and several values of incidence

Description

Calculates negative log-likelihood (deviance) for one antigen:isotype pair and several values of incidence (lambda).

Usage

f_dev(lambda, csdata, lnpars, cond)

Arguments

Details

Vectorized version of f_dev0(); interface with C lib serocalc.so

Value

a numeric vector of negative log-likelihoods, corresponding to the elements of input lambda

Examples

library(dplyr)
library(tibble)

# load in longitudinal parameters
curve_params <-
  typhoid_curves_nostrat_100 %>%
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

# load in pop data
xs_data <-
  sees_pop_data_pk_100

#Load noise params
noise_params <- 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)

cur_antibody = "HlyE_IgA"

cur_data =
  xs_data %>%
  dplyr::filter(
   .data$catchment == "aku",
   .data$antigen_iso == cur_antibody) %>%
  dplyr::slice_head(n = 100)

cur_curve_params =
  curve_params %>%
  dplyr::filter(.data$antigen_iso == cur_antibody) %>%
  dplyr::slice_head(n = 100)

cur_noise_params =
  noise_params %>%
  dplyr::filter(.data$antigen_iso == cur_antibody)

if(!is.element('d', names(cur_curve_params)))
{
  cur_curve_params =
    cur_curve_params %>%
    dplyr::mutate(
      alpha = .data$alpha * 365.25,
      d = .data$r - 1)
}

lambdas = seq(.1, .2, by = .01)
f_dev(
    lambda = lambdas,
    csdata = cur_data,
    lnpars = cur_curve_params,
    cond = cur_noise_params
  )

Calculate negative log-likelihood (deviance) for one antigen:isotype pair and one incidence rate

Description

Calculate negative log-likelihood (deviance) for one antigen:isotype pair and one incidence rate

Usage

f_dev0(lambda, csdata, lnpars, cond)

Arguments

Details

interface with C lib serocalc.so

Value

a numeric() negative log-likelihood, corresponding to input lambda

Examples

library(dplyr)
library(tibble)

# load in longitudinal parameters
curve_params <-
  typhoid_curves_nostrat_100 %>%
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

# load in pop data
xs_data <-
  sees_pop_data_pk_100

#Load noise params
noise_params <- 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)

cur_antibody = "HlyE_IgA"

cur_data <-
  xs_data %>%
  dplyr::filter(
   .data$antigen_iso == cur_antibody) %>%
  dplyr::slice_head(n = 100)

cur_curve_params <-
  curve_params %>%
  dplyr::filter(.data$antigen_iso == cur_antibody) %>%
  dplyr::slice_head(n = 100)

cur_noise_params <-
  noise_params %>%
  dplyr::filter(.data$antigen_iso == cur_antibody)

if (!is.element('d', names(cur_curve_params)))
{
  cur_curve_params <-
    cur_curve_params %>%
    dplyr::mutate(
      alpha = .data$alpha * 365.25,
      d = .data$r - 1)
}

lambda = 0.1
f_dev0(
    lambda = lambda,
    csdata = cur_data,
    lnpars = cur_curve_params,
    cond = cur_noise_params
  )

Calculate negative log-likelihood (deviance)

Description

fdev() was renamed to f_dev() to create a more consistent API.

Usage

fdev(lambda, csdata, lnpars, cond)

Extract biomarker levels

Description

Extract biomarker levels

Usage

get_biomarker_levels(object, ...)

Arguments

Value

the biomarker levels in object

Examples

sees_pop_data_100 |> get_biomarker_levels()

Extract biomarker names column

Description

Extract biomarker names column

Usage

get_biomarker_names(object, ...)

Arguments

Value

a character or factor vector of biomarker names

Examples

sees_pop_data_100 |> get_biomarker_names() |> head()

Get biomarker variable name

Description

Get biomarker variable name

Usage

get_biomarker_names_var(object, ...)

Arguments

Value

a character string identifying the biomarker names column in object

Examples

sees_pop_data_100 |> get_biomarker_names_var()

Get antibody measurement values

Description

Get antibody measurement values

Usage

get_values(object, ...)

Arguments

Value

a numeric vector of antibody measurement values

Examples

sees_pop_data_100 |> get_values()

Extract antibody measurement values

Description

Extract antibody measurement values

Usage

get_values_var(object, ...)

Arguments

Value

the name of the column in object specified as containing antibody abundance measurements

Examples

sees_pop_data_100 |> get_values_var()

Graph estimated antibody decay curves

Description

Graph estimated antibody decay curves

Usage

graph.curve.params(
  object,
  antigen_isos = unique(object$antigen_iso),
  verbose = FALSE,
  quantiles = c(0.1, 0.5, 0.9),
  alpha_samples = 0.3,
  chain_color = TRUE,
  log_x = FALSE,
  log_y = TRUE,
  n_curves = 100,
  iters_to_graph = head(unique(object$iter), n_curves),
  ...
)

Arguments

Details

n_curves and iters_to_graph

In most cases, object 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 iters_to_graph argument, then n_curves has no effect.

Value

a ggplot2::ggplot() object showing the antibody dynamic kinetics of selected antigen/isotype combinations, with optional posterior distribution quantile curves.

Examples

# Load example dataset
curve <- typhoid_curves_nostrat_100 |>
  dplyr::filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

# Plot quantiles without showing all curves
plot1 <- graph.curve.params(curve, n_curves = 0)
print(plot1)

# Plot with additional quantiles and show all curves
plot2 <- graph.curve.params(
  curve,
  n_curves = Inf,
  quantiles = c(0.1, 0.5, 0.9)
)
print(plot2)

# Plot with MCMC chains in black
plot3 <- graph.curve.params(
  curve,
  n_curves = Inf,
  quantiles = c(0.1, 0.5, 0.9),
  chain_color = FALSE
)
print(plot3)

Graph log-likelihood of data

Description

Graph log-likelihood of data

Usage

graph_loglik(
  pop_data,
  curve_params,
  noise_params,
  antigen_isos = pop_data %>% get_biomarker_levels(),
  x = 10^seq(-3, 0, by = 0.1),
  highlight_points = NULL,
  highlight_point_names = "highlight_points",
  log_x = FALSE,
  previous_plot = NULL,
  curve_label = paste(antigen_isos, collapse = " + "),
  ...
)

Arguments

Value

a ggplot2::ggplot()

Examples

library(dplyr)
library(tibble)

# Load cross-sectional data
xs_data <-
  sees_pop_data_pk_100

# Load curve parameters and subset for the purposes of this example
curve <-
  typhoid_curves_nostrat_100 %>%
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

# Load noise parameters
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)

# Graph the log likelihood
lik_HlyE_IgA <- # nolint: object_name_linter
  graph_loglik(
    pop_data = xs_data,
    curve_params = curve,
    noise_params = cond,
    antigen_isos = "HlyE_IgA",
    log_x = TRUE
)

lik_HlyE_IgA # nolint: object_name_linter

graph antibody decay curves by antigen isotype

Description

graph antibody decay curves by antigen isotype

Usage

graph_seroresponse_model_1(
  object,
  antigen_isos = unique(object$antigen_iso),
  ncol = min(3, length(antigen_isos)),
  ...
)

Arguments

Details

iters_to_graph

If you directly specify iters_to_graph when calling this function, the row numbers are enumerated separately for each antigen isotype; in other words, for the purposes of this argument, row numbers start over at 1 for each antigen isotype. There is currently no way to specify different row numbers for different antigen isotypes; if you want to do that, you will could call plot_curve_params_one_ab() directly for each antigen isotype and combine the resulting panels yourself. Or you could subset curve_params manually, before passing it to this function, and set the n_curves argument to Inf.

Value

a ggplot2::ggplot() object

Examples

library(dplyr)
library(ggplot2)
library(magrittr)

curve <-
  serocalculator_example("example_curve_params.csv") |>
  read.csv() |>
  as_sr_params() |>
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG")) |>
  graph_seroresponse_model_1()

curve

Get person IDs

Description

Get person IDs

Usage

ids(object, ...)

Arguments

Value

a character vector

Examples

ids(sees_pop_data_pk_100)

Get person ID variable name

Description

Get person ID variable name

Usage

ids_varname(object, ...)

Arguments

Value

a character string containing the person ID column, as recorded in the metadata of object

Examples

ids_varname(sees_pop_data_pk_100)

extract a row from longitudinal parameter set

Description

take a random sample from longitudinal parameter set given age at infection, for a list of antibodies

Usage

ldpar(age, antigen_isos, nmc, npar, ...)

Arguments

Value

an array of 7 parameters and initial conditions: c(y0,b0,mu0,mu1,c1,alpha,shape)

Calculate log-likelihood

Description

Calculates the log-likelihood of a set of cross-sectional antibody response data, for a given incidence rate (lambda) value.

llik() was renamed to log_likelihood() to create a more consistent API.

Usage

llik(...)

Load antibody decay curve parameter samples

Description

load_curve_params() was renamed to load_sr_params() to create a more consistent API.

Usage

load_curve_params(...)

Load noise parameters

Description

Load noise parameters

Usage

load_noise_params(file_path, antigen_isos = NULL)

Arguments

Value

a noise object (a tibble::tbl_df with extra attribute antigen_isos)

Examples

noise <- load_noise_params(serocalculator_example("example_noise_params.rds"))
print(noise)

Load a cross-sectional antibody survey data set

Description

Load a cross-sectional antibody survey data set

Usage

load_pop_data(file_path, ...)

Arguments

Value

a pop_data object (a tibble::tbl_df with extra attributes)

Examples

xs_data <- load_pop_data(serocalculator_example("example_pop_data.rds"))

print(xs_data)

Load longitudinal seroresponse parameter samples

Description

Load longitudinal seroresponse parameter samples

Usage

load_sr_params(file_path, antigen_isos = NULL)

Arguments

Value

a curve_params object (a tibble::tbl_df with extra attribute antigen_isos)

Examples

curve <- load_sr_params(serocalculator_example("example_curve_params.rds"))

print(curve)

Calculate log-likelihood

Description

Calculates the log-likelihood of a set of cross-sectional antibody response data, for a given incidence rate (lambda) value.

Usage

log_likelihood(
  lambda,
  pop_data,
  curve_params,
  noise_params,
  antigen_isos = get_biomarker_levels(pop_data),
  verbose = FALSE,
  ...
)

Arguments

Value

the log-likelihood of the data with the current parameter values

Examples

library(dplyr)
library(tibble)

# Load cross-sectional data
xs_data <-
  sees_pop_data_pk_100

# Load curve parameters and subset for the purposes of this example
curve <-
  typhoid_curves_nostrat_100 %>%
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

# 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(
  pop_data = xs_data,
  curve_params = curve,
  noise_params = cond,
  antigen_isos = c("HlyE_IgG", "HlyE_IgA"),
  lambda = 0.1
) %>% print()

generate random sample from baseline distribution

Description

generate random sample from baseline distribution

Usage

mk_baseline(kab, n = 1, blims, ...)

Arguments

Value

a numeric() vector

Graph an antibody decay curve model

Description

Graph an antibody decay curve model

Usage

plot_curve_params_one_ab(
  object,
  verbose = FALSE,
  alpha = 0.4,
  n_curves = 100,
  n_points = 1000,
  log_x = FALSE,
  log_y = TRUE,
  iters_to_graph = seq_len(min(n_curves, nrow(object))),
  xlim = c(10^-1, 10^3.1),
  ...
)

Arguments

Details

n_curves and iters_to_graph

In most cases, object 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 iters_to_graph argument, then n_curves has no effect.

Value

a ggplot2::ggplot() object

Examples

library(dplyr) # loads the `%>%` operator and `dplyr::filter()`

curve <-
  typhoid_curves_nostrat_100 %>%
  filter(antigen_iso == ("HlyE_IgG")) %>%
  serocalculator:::plot_curve_params_one_ab()

  curve

Print Method for seroincidence Class

Description

print() function for seroincidence objects from est_seroincidence()

Usage

## S3 method for class 'seroincidence'
print(x, ...)

Arguments

Value

an invisible copy of input parameter x

Examples

library(dplyr)

xs_data <-
  sees_pop_data_pk_100

curve <-
  typhoid_curves_nostrat_100 %>%
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

noise <-
  example_noise_params_pk

est1 <- est_seroincidence(
  pop_data = xs_data,
  sr_params = curve,
  noise_params = noise,
  antigen_isos = c("HlyE_IgG", "HlyE_IgA"),
)
print(est1)

Print Method for seroincidence.by Object

Description

Custom print() function for seroincidence.by objects (from est_seroincidence_by())

Usage

## S3 method for class 'seroincidence.by'
print(x, ...)

Arguments

Value

an invisible copy of input parameter x

Examples

library(dplyr)

xs_data <-
  sees_pop_data_pk_100

curve <-
  typhoid_curves_nostrat_100 |>
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

noise <-
  example_noise_params_pk

# estimate seroincidence
est2 <- est_seroincidence_by(
  strata = c("catchment"),
  pop_data = xs_data,
  sr_params = curve,
  noise_params = noise,
  antigen_isos = c("HlyE_IgG", "HlyE_IgA"),
  # num_cores = 8 # Allow for parallel processing to decrease run time
)

# calculate summary statistics for the seroincidence object
print(est2)

Print Method for Seroincidence Summary Object

Description

Custom print() function for "summary.seroincidence.by" objects (constructed by summary.seroincidence.by())

Usage

## S3 method for class 'summary.seroincidence.by'
print(x, ...)

Arguments

Value

an invisible copy of input parameter x

Examples

library(dplyr)

xs_data <-
  sees_pop_data_pk_100

curve <-
  typhoid_curves_nostrat_100 |>
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

noise <-
  example_noise_params_pk

# estimate seroincidence
est2 <- est_seroincidence_by(
  strata = c("catchment"),
  pop_data = xs_data,
  sr_params = curve,
  noise_params = noise,
  antigen_isos = c("HlyE_IgG", "HlyE_IgA"),
  # num_cores = 8 # Allow for parallel processing to decrease run time
)

# calculate summary statistics for the seroincidence object
est2_summary <- summary(est2)
print(est2_summary)

Objects exported from other packages

Description

These objects are imported from other packages. Follow the links below to see their documentation.

ggplot2

autoplot

Small example cross-sectional data set

Description

A subset of data from the SEES data, for examples and testing.

Usage

sees_pop_data_100

Format

sees_pop_data_pk_100

A pop_data object (from as_pop_data()) with 200 rows and 8 columns:

id

Observation ID

Country

Country where the participant was living

cluster

survey sampling cluster

catchment

survey catchment area

age

participant's age when sampled, in years

antigen_iso

which antigen and isotype are being measured (data is in long format)

value

concentration of antigen isotype, in ELISA units

Source

https://osf.io/n6cp3

Small example cross-sectional data set

Description

A subset of data from the SEES data, for examples and testing, data from Pakistan only.

Usage

sees_pop_data_pk_100

Format

sees_pop_data_pk_100

A pop_data object (from as_pop_data()) with 200 rows and 8 columns:

id

Observation ID

Country

Country where the participant was living

cluster

survey sampling cluster

catchment

survey catchment area

age

participant's age when sampled, in years

antigen_iso

which antigen and isotype are being measured (data is in long format)

value

concentration of antigen isotype, in ELISA units

Source

https://osf.io/n6cp3

Small example cross-sectional data set

Description

A subset of data from the SEES data, for examples and testing, data from Pakistan only, variable names not normalized by as_pop_data().

Usage

sees_pop_data_pk_100_old_names

Format

sees_pop_data_pk_100_old_names

A pop_data object (from as_pop_data()) with 200 rows and 8 columns:

index_id

Observation ID

Country

Country where the participant was living

cluster

survey sampling cluster

catchment

survey catchment area

Age

participant's age when sampled, in years

antigen_iso

which antigen and isotype are being measured (data is in long format)

result

concentration of antigen isotype, in ELISA units

Source

https://osf.io/n6cp3

Example "seroincidence.by" object

Description

Typhoid seroconversion rate estimates by country and age category from the SEES study.

Usage

sees_typhoid_ests_strat

Format

An object of class seroincidence.by (inherits from list) of length 9.

Source

serocalculator/data-raw/sees_typhoid_ests_strat.R

Get path to an example file

Description

The serocalculator package comes bundled with a number of sample files in its inst/extdata directory. This serocalculator_example() function make those sample files easy to access.

Usage

serocalculator_example(file = NULL)

Arguments

Details

Adapted from readr::readr_example() following the guidance in https://r-pkgs.org/data.html#sec-data-example-path-helper.

Value

a character string providing the path to the file specified by file, or a vector or available files if file = NULL.

Examples

serocalculator_example()
serocalculator_example("example_pop_data.csv")

Specify person ID column

Description

Sets the id_var metadata attribute of object

Usage

set_id_var(object, id = "index_id", standardize = TRUE, ...)

Arguments

Value

a modified version of object

Examples

serocalculator_example("example_pop_data.rds") |>
  readr::read_rds() |>
    set_id_var(id = "index_id") |>
    attr("id_var")

Simulate a cross-sectional serosurvey with noise

Description

sim.cs() was renamed to sim_pop_data() to create a more consistent API.

Usage

sim.cs(n.smpl, age.rng, n.mc, renew.params, add.noise, ...)

Simulate multiple data sets

Description

sim.cs.multi() was renamed to sim_pop_data_multi() to create a more consistent API.

Usage

sim.cs.multi(...)

Simulate a cross-sectional serosurvey with noise

Description

Makes a cross-sectional data set (age, y(t) set) and adds noise, if desired.

Usage

sim_pop_data(
  lambda = 0.1,
  n_samples = 100,
  age_range = c(0, 20),
  age_fixed = NA,
  antigen_isos = intersect(get_biomarker_levels(curve_params), rownames(noise_limits)),
  n_mcmc_samples = 0,
  renew_params = FALSE,
  add_noise = FALSE,
  curve_params,
  noise_limits,
  format = "wide",
  verbose = FALSE,
  ...
)

Arguments

Value

a tibble::tbl_df containing simulated cross-sectional serosurvey data, with columns:

• age: age (in days)

• one column for each element in the antigen_iso input argument

Examples

# Load curve parameters
dmcmc <- typhoid_curves_nostrat_100

# Specify the antibody-isotype responses to include in analyses
antibodies <- c("HlyE_IgA", "HlyE_IgG")

# Set seed to reproduce results
set.seed(54321)

# Simulated incidence rate per person-year
lambda <- 0.2
# Range covered in simulations
lifespan <- c(0, 10)
# Cross-sectional sample size
nrep <- 100

# Biologic noise distribution
dlims <- rbind(
  "HlyE_IgA" = c(min = 0, max = 0.5),
  "HlyE_IgG" = c(min = 0, max = 0.5)
)

# Generate cross-sectional data
csdata <- sim_pop_data(
  curve_params = dmcmc,
  lambda = lambda,
  n_samples = nrep,
  age_range = lifespan,
  antigen_isos = antibodies,
  n_mcmc_samples = 0,
  renew_params = TRUE,
  add_noise = TRUE,
  noise_limits = dlims,
  format = "long"
)

Simulate multiple data sets

Description

Simulate multiple data sets

Usage

sim_pop_data_multi(
  nclus = 10,
  sample_sizes = 100,
  lambdas = c(0.05, 0.1, 0.15, 0.2, 0.3),
  num_cores = max(1, parallel::detectCores() - 1),
  rng_seed = 1234,
  verbose = FALSE,
  ...
)

Arguments

Value

a tibble::tibble()

Examples

# Load curve parameters
dmcmc <- typhoid_curves_nostrat_100

# Specify the antibody-isotype responses to include in analyses
antibodies <- c("HlyE_IgA", "HlyE_IgG")

# Set seed to reproduce results
set.seed(54321)

# Simulated incidence rate per person-year
lambdas = c(.05, .1, .15, .2, .3)

# Range covered in simulations
lifespan <- c(0, 10);

# Cross-sectional sample size
nrep <- 100

# Biologic noise distribution
dlims <- rbind(
  "HlyE_IgA" = c(min = 0, max = 0.5),
  "HlyE_IgG" = c(min = 0, max = 0.5)
)

sim_data <- sim_pop_data_multi(
  curve_params = dmcmc,
  lambdas = lambdas,
  sample_sizes = nrep,
  age_range = lifespan,
  antigen_isos = antibodies,
  n_mcmc_samples = 0,
  renew_params = TRUE,
  add_noise = TRUE,
  noise_limits = dlims,
  format = "long",
  nclus = 10)

sim_data

collect cross-sectional data

Description

collect cross-sectional data

Usage

simcs.tinf(
  lambda,
  n_samples,
  age_range,
  age_fixed = NA,
  antigen_isos,
  n_mcmc_samples = 0,
  renew_params = FALSE,
  ...
)

Arguments

Value

an array() with dimensions n_samples, length(antigen_isos) + 1, where rows are observations and columns are age and biomarkers y(t)

simulate antibody kinetics of y over a time interval

Description

simulate antibody kinetics of y over a time interval

Usage

simresp.tinf(
  lambda,
  t_end,
  age_fixed,
  antigen_isos,
  n_mcmc_samples = 0,
  renew_params,
  predpar,
  ...
)

Arguments

Value

a list with:

• t = times (in days, birth at day 0),

• b = bacteria level, for each antibody signal (not used; probably meaningless),

• y = antibody level, for each antibody signal

• smp = whether an infection involves a big jump or a small jump

• t.inf = times when infections have occurred.

Barplot method for summary.seroincidence.by objects

Description

Barplot method for summary.seroincidence.by objects

Usage

strat_ests_barplot(
  object,
  yvar,
  color_var = NULL,
  alpha = 0.7,
  CIs = FALSE,
  title = NULL,
  xlab = "Seroconversion rate per 1000 person-years",
  ylab = yvar,
  fill_lab = NULL,
  color_palette = NULL,
  ...
)

Arguments

Value

a ggplot2::ggplot() object.

Scatterplot method for summary.seroincidence.by objects

Description

Scatterplot method for summary.seroincidence.by objects

Usage

strat_ests_scatterplot(
  object,
  xvar = strata(object)[1],
  alpha = 0.7,
  shape = 1,
  dodge_width = 0.001,
  CIs = FALSE,
  color_var = "nlm.convergence.code",
  group_var = NULL,
  ...
)

Arguments

Value

a ggplot2::ggplot() object

Examples

library(dplyr)
library(ggplot2)

xs_data <-
  sees_pop_data_pk_100

curve <-
  typhoid_curves_nostrat_100 |>
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

noise <-
  example_noise_params_pk

est2 <- est_seroincidence_by(
  strata = c("catchment", "ageCat"),
  pop_data = xs_data,
  sr_params = curve,
  noise_params = noise,
  curve_strata_varnames = NULL,
  noise_strata_varnames = NULL,
  antigen_isos = c("HlyE_IgG", "HlyE_IgA"),
  num_cores = 2 # Allow for parallel processing to decrease run time
)

est2sum <- summary(est2)

strat_ests_scatterplot(est2sum,
                       xvar = "ageCat",
                       color_var = "catchment",
                       CIs = TRUE,
                       group_var = "catchment")

Extract Strata metadata from an object

Description

Generic method for extracting strata metadata from objects. See strata.default()

Usage

strata(x)

Arguments

Value

the strata metadata of x

Extract the Strata attribute from an object, if present

Description

Extract the Strata attribute from an object, if present

Usage

## Default S3 method:
strata(x)

Arguments

Value

• a tibble::tibble() with strata in rows, or

• NULL if x does not have a "strata" attribute

Split data by stratum

Description

Split biomarker data, decay curve parameters, and noise parameters to prepare for stratified incidence estimation.

Usage

stratify_data(
  data,
  curve_params,
  noise_params,
  strata_varnames = "",
  curve_strata_varnames = NULL,
  noise_strata_varnames = NULL,
  antigen_isos = get_biomarker_levels(data)
)

Arguments

Value

a "biomarker_data_and_params.list" object (a list with extra attributes "strata", "antigen_isos", etc)

Examples

## Not run: 
library(dplyr)

xs_data <-
  sees_pop_data_pk_100

curve <-
  typhoid_curves_nostrat_100 |>
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

noise <-
  example_noise_params_pk

stratified_data =
  stratify_data(
   data = xs_data,
   curve_params = curve,
   noise_params = noise,
   strata_varnames = "catchment",
   curve_strata_varnames = NULL,
   noise_strata_varnames = NULL
   )

## End(Not run)

Summarize cross-sectional antibody survey data

Description

summary() method for pop_data objects

Usage

## S3 method for class 'pop_data'
summary(object, strata = NULL, ...)

## S3 method for class 'summary.pop_data'
print(x, ...)

Arguments

Value

a summary.pop_data object, which is a list containing two summary tables:

• age_summary summarizing age

• ab_summary summarizing value, stratified by antigen_iso

Examples

library(dplyr)

xs_data <-
  sees_pop_data_pk_100
summary(xs_data, strata = "catchment")

Summarizing fitted seroincidence models

Description

This function is a summary() method for seroincidence objects.

Usage

## S3 method for class 'seroincidence'
summary(object, coverage = 0.95, verbose = TRUE, ...)

Arguments

Value

a tibble::tibble() containing the following:

• est.start: the starting guess for incidence rate

• ageCat: the age category we are analyzing

• incidence.rate: the estimated incidence rate, per person year

• CI.lwr: lower limit of confidence interval for incidence rate

• CI.upr: upper limit of confidence interval for incidence rate

• coverage: coverage probability

• log.lik: log-likelihood of the data used in the call to est_seroincidence(), evaluated at the maximum-likelihood estimate of lambda (i.e., at incidence.rate)

• iterations: the number of iterations used

• antigen_isos: a list of antigen isotypes used in the analysis

• nlm.convergence.code: information about convergence of the likelihood maximization procedure performed by nlm() (see "Value" section of stats::nlm(), component code); codes 3-5 indicate issues:

  • 1: relative gradient is close to zero, current iterate is probably solution.

  • 2: successive iterates within tolerance, current iterate is probably solution.

  • 3: Last global step failed to locate a point lower than x. Either x is an approximate local minimum of the function, the function is too non-linear for this algorithm, or stepmin in est_seroincidence() (a.k.a., steptol in stats::nlm()) is too large.

  • 4: iteration limit exceeded; increase iterlim.

  • 5: maximum step size stepmax exceeded five consecutive times. Either the function is unbounded below, becomes asymptotic to a finite value from above in some direction, or stepmax is too small.

Examples

library(dplyr)

xs_data <-
  sees_pop_data_pk_100

curve <-
  typhoid_curves_nostrat_100 |>
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

noise <-
  example_noise_params_pk

est1 <- est_seroincidence(
  pop_data = xs_data,
  sr_params = curve,
  noise_params = noise,
  antigen_isos = c("HlyE_IgG", "HlyE_IgA")
)

summary(est1)

Summary Method for "seroincidence.by" Objects

Description

Calculate seroincidence from output of the seroincidence calculator est_seroincidence_by().

Usage

## S3 method for class 'seroincidence.by'
summary(
  object,
  confidence_level = 0.95,
  show_deviance = TRUE,
  show_convergence = TRUE,
  verbose = FALSE,
  ...
)

Arguments

Value

A summary.seroincidence.by object, which is a tibble::tibble, with the following columns:

• incidence.rate maximum likelihood estimate of lambda (seroincidence)

• CI.lwr lower confidence bound for lambda

• CI.upr upper confidence bound for lambda

• Deviance (included if show_deviance = TRUE) Negative log likelihood (NLL) at estimated (maximum likelihood) lambda)

• nlm.convergence.code (included if show_convergence = TRUE) Convergence information returned by stats::nlm()

The object also has the following metadata (accessible through base::attr()):

• antigen_isos Character vector with names of input antigen isotypes used in est_seroincidence_by()

• Strata Character with names of strata used in est_seroincidence_by()

Examples

library(dplyr)

xs_data <-
  sees_pop_data_pk_100

curve <-
  typhoid_curves_nostrat_100 |>
  filter(antigen_iso %in% c("HlyE_IgA", "HlyE_IgG"))

noise <-
  example_noise_params_pk

# estimate seroincidence
est2 <- est_seroincidence_by(
  strata = c("catchment"),
  pop_data = xs_data,
  sr_params = curve,
  noise_params = noise,
  antigen_isos = c("HlyE_IgG", "HlyE_IgA"),
  # num_cores = 8 # Allow for parallel processing to decrease run time
)

# calculate summary statistics for the seroincidence object
summary(est2)

Small example of antibody response curve parameters for typhoid

Description

A subset of data from the SEES study, for examples and testing.

Usage

typhoid_curves_nostrat_100

Format

typhoid_curves_nostrat_100

A curve_params object (from as_sr_params()) with 500 rows and 7 columns:

antigen_iso

which antigen and isotype are being measured (data is in long format)

iter

MCMC iteration

y0

Antibody concentration at t = 0 (start of active infection)

y1

Antibody concentration at t = t1 (end of active infection)

t1

Duration of active infection

alpha

Antibody decay rate coefficient

r

Antibody decay rate exponent parameter

Source

https://osf.io/rtw5k

Warn about missing stratifying variables in a dataset

Description

Warn about missing stratifying variables in a dataset

Usage

warn_missing_strata(data, strata, dataname)

Arguments

Value

a character() vector of the subset of stratifying variables that are present in pop_data

Examples

## Not run: 
expected_strata <- data.frame(Species = "banana", type = "orchid")

warn_missing_strata(iris, expected_strata, dataname = "iris")

## End(Not run)