truePrevPools: Estimate true prevalence from pooled samples

Description

Bayesian estimation of true prevalence from apparent prevalence obtained by testing pooled samples.

Usage

truePrevPools(x, n, SE = 1, SP = 1, prior = c(1, 1),
              nchains = 2, burnin = 10000, update = 10000,
              verbose = FALSE)

Arguments

The vector of indicator variables, indicating whether a pool was positive ("1") or negative ("0")

The vector of pool sizes

SE, SP

The prior distribution for sensitivity (SE) and specificity (SP); see 'Details' below for specification of these distributions

prior

The parameters of the prior Beta distribution for true prevalence; defaults to c(1, 1)

nchains

The number of chains used in the estimation process; nchains must be \(\ge 2\)

burnin

The number of discarded model iterations; defaults to 10,000

update

The number of withheld model iterations; defaults to 10,000

verbose

Logical flag, indicating if JAGS process output should be printed to the R console; defaults to FALSE

Value

An object of class prev.

Details

truePrevPools calls on JAGS/rjags to estimate the true prevalence from the apparent prevalence in a Bayesian framework. The default model, in BUGS language, is given below. To see the actual fitted model, see the model slot of the prev-object.

  model {
    for (i in 1:N) {
      x[i] ~ dbern(AP[i])
      AP[i] <- SEpool[i] * (1 - pow(1 - TP, n[i])) + (1 - SPpool[i]) * pow(1 - TP, n[i])
      SEpool[i] <- 1 - (pow(1 - SE, n[i] * TP) * pow(SP, n[i] * (1 - TP)))
      SPpool[i] <- pow(SP, n[i])
    }
  # SE ~ user-defined (see below)
  # SP ~ user-defined (see below)
  TP ~ dbeta(prior[1], prior[2])
  }

The test sensitivity (SE) and specificity (SP) can be specified by the user, independently, as one of "fixed", "uniform", "beta", "pert", or "beta-expert", with "fixed" as the default. Note that SE and SP must correspond to the test characteristics for testing individual samples; truePrevPools will calculate SEpool and SPpool, the sensitivity and specificitiy for testing pooled samples, based on Boelaert et al. (2000).

Distribution parameters can be specified in a named list() as follows:

Fixed: list(dist = "fixed", par)
Uniform: list(dist = "uniform", min, max)
Beta: list(dist = "beta", alpha, beta)
PERT: list(dist = "pert", method, a, m, b, k) 'method' must be "classic" or "vose"; 'a' denotes the pessimistic (minimum) estimate, 'm' the most likely estimate, and 'b' the optimistic (maximum) estimate; 'k' denotes the scale parameter. See betaPERT for more information on Beta-PERT parametrization.
Beta-Expert: list(dist = "beta-expert", mode, mean, lower, upper, p) 'mode' denotes the most likely estimate, 'mean' the mean estimate; 'lower' denotes the lower bound, 'upper' the upper bound; 'p' denotes the confidence level of the expert. Only mode or mean should be specified; lower and upper can be specified together or alone. See betaExpert for more information on Beta-Expert parameterization.

For Uniform, Beta and Beta-PERT distributions, BUGS-style short-hand notation is also allowed:

Uniform: ~dunif(min, max)
Beta: ~dbeta(alpha, beta)
Beta-PERT: ~dpert(min, mode, max)

References

Speybroeck N, Williams CJ, Lafia KB, Devleesschauwer B, Berkvens D (2012) Estimating the prevalence of infections in vector populations using pools of samples. Med Vet Entomol 26:361-371
Boelaert F, Walravens K, Biront P, Vermeersch JP, Berkvens D, Godfroid J (2000) Prevalence of paratuberculosis (Johne's disease) in the Belgian cattle population. Vet Microbiol 77:269-281

Examples

Run this code

# NOT RUN {
## Sandflies in Aurabani, Nepal, 2007
pool_results <- c(0, 0, 0, 0, 0, 0, 0, 0, 1, 0)
pool_sizes <- c(2, 1, 6, 10, 1, 7, 1, 4, 1, 3)

## Sensitivity ranges uniformly between 60% and 95%
## Specificity is considered to be 100%

#> BUGS-style:
truePrevPools(x = pool_results, n = pool_sizes,
              SE = ~dunif(0.60, 0.95), SP = 1)
			  
#> list-style:
SE <- list(dist = "uniform", min = 0.60, max = 0.95)
truePrevPools(x = pool_results, n = pool_sizes,
              SE = SE, SP = 1)
# }

Run the code above in your browser using DataLab