In ku-awdc/eggSim: Simulating Datasets Representing Egg Reduction Rate Data

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

survey_sim function

TODO: tidy up and expand

Installation

The eggSim package is hosted via a public repository at (https://github.com/ku-awdc/eggSim), so you can install it directly from there using remotes::install_github(). However, this requires rebuilding of vignettes and always requires C++ compilers. We therefore recommend that you install from our drat repository using the following command:

install.packages('eggSim', repos=c(CRAN="https://cran.rstudio.com/", 
                                   "ku-awdc"="https://ku-awdc.github.io/drat/"))

This will install a pre-compiled stable/release version where this is available for your system and R version, but may still ask for installation from source if there is a newer source version compared to the pre-built binary that is available for your platform. To be able to install these you need to do one of the following:

On macOS run the following command from the terminal: sudo xcode-select --install
On Windows download and install Rtools (make sure to pick the version that matches your R installation) from: https://cran.r-project.org/bin/windows/Rtools/

Running the code shown in this vignette requires at least version 0.9.4-6:

stopifnot(packageVersion("eggSim") >= "0.9.4-6")

Introduction

The eggSim pacakge is intended to facilitate simulation-based assessment of survey design strategies for testing anthelmintic efficacy via faecal egg count reduction (FECRT) tests. The software is currently a work in progress, and comments/suggestions/contributions are very welcome to the package maintainers!

The package is designed to work in a tidyverse-like way, so loading the tidyverse meta-package along with the eggSim package is a good idea:

library("eggSim")
library("tidyverse")

Scenarios

First define the scenarios to test. This should be a data frame with columns parasite, scenario, mean_epg, true_efficacy and cutoff. It is most convenient to use the provided helper function to generate a data frame in the correct format:

scenarios <- survey_scenario(parasite = c("ascaris","trichuris","hookworm"))
scenarios

If needed, you can then modify this e.g.:

scenarios$mean_epg[1] <- 20
scenarios

Then you need to define the parameters, either as a data frame for a single parameter set, or a list of data frames for multiple parameter sets. Each data frame can either have 1 row, or multiple rows corresponding to iterations where uncertainty in parameter estimates will be integrated over by the simulation. Again, it is convinient to use the provided helper function:

par_1 <- survey_parameters(design = c("NS_11"), parasite = c("ascaris"), method = c("kk"))
par_1

par_all <- survey_parameters(design = c("SS_11","SS_12","NS_11","NS_12","SSR_11","SSR_12"), parasite = c("ascaris","trichuris","hookworm"), method = c("kk","miniflotac","fecpak"))
length(par_all)
par_all[[1]]

Again, the values can be changed e.g.:

par_all[[1]]$reduction_cv <- 1

Or if you want to integrate over uncertainty in one or more parameter estimate, you can alter the data frame so that it has as many rows as iterations, with each row of the data frame reflecting uncertainty in one or more of the parameter estimates (which is most likely to apply to the cv estimates). For example:

par_all[[2]] <- par_all[[2]] |> slice(rep(1, 1e3)) |> mutate(reduction_cv = rgamma(n(), 1, 1))
nrow(par_all[[2]])

Then you can run the simulation, specifying a vector of number of individuals, number of iterations, and optionally the number of parallel cores to use. Here I select only the first 5 parameter sets to speed up processing of the vignette:

results <- survey_sim(n_individ = seq(100,1000,by=10), scenario = scenarios, parameters = par_all[1:5], iterations = 1e3, cl=2L, output="summarised")

Note that the default option for the output is to pre-summarise (in order to reduce memory usage). If you prefer, then you can ask for either full output or extended output, which return results at the level of iteration:

results_full <- survey_sim(n_individ = seq(100,1000,by=10), scenario = scenarios, parameters = par_1, iterations = 1e3, cl=2L, output="full")

Here are some suggested plots:

pp <- "ascaris"
ggplot(results |> filter(parasite==pp), aes(x=n_individ, y=(n_failure/n_total), col=design, group=design)) +
  geom_line() +
  facet_grid(mean_epg ~ method) +
  labs(col=pp)

ggplot(results |> filter(parasite==pp), aes(x=n_individ, y=n_below_cutoffs/n_total, col=design, group=design)) +
    geom_line() +
    facet_grid(mean_epg ~ method) +
    labs(col=pp)

ggplot(results |> filter(parasite==pp), aes(x=-cost_mean, y=n_below_cutoffs/n_total, col=design, group=design)) +
  geom_line() +
  facet_grid(mean_epg ~ method) +
  labs(col=pp)

ggplot(results |> filter(parasite==pp), aes(x=-cost_mean, y=n_below_cutoffs/n_total, col=design, group=design)) +
  geom_line() +
  facet_grid(mean_epg ~ method) +
  labs(col=pp) +
  coord_cartesian(ylim=c(0.5, NA))

Contributing

We would highly value any contributions (either comments, suggestions or code patches) on this package! You can contribute by either working with the GitHub repository or by emailing the package maintainers using the email addresses available from: