Unit tests (#71)

* Updated plot_sensitivity

* Added example for plot_sensitivity

* Added unit test for GeoTox S3 object

R/GeoTox.R

@@ -5,6 +5,46 @@
 #' GeoTox analysis.
 #'
 #' @export
+#' 
+#' @examples
+#' # Use a subset of the package data for demonstration purposes
+#' set.seed(2357)
+#' n <- 10 # Population size
+#' m <- 5 # Number of regions
+#' idx <- if (m < 100) sample(1:100, m) else 1:100
+#' 
+#' geoTox <- GeoTox() |> 
+#'   # Set region and group boundaries (for plotting)
+#'   set_boundaries(region = geo_tox_data$boundaries$county,
+#'                  group  = geo_tox_data$boundaries$state) |> 
+#'   # Simulate populations for each region
+#'   simulate_population(age           = split(geo_tox_data$age, ~FIPS)[idx],
+#'                       obesity       = geo_tox_data$obesity[idx, ],
+#'                       exposure      = split(geo_tox_data$exposure, ~FIPS)[idx],
+#'                       simulated_css = geo_tox_data$simulated_css,
+#'                       n             = n) |> 
+#'   # Estimated Hill parameters
+#'   set_hill_params(geo_tox_data$dose_response |>
+#'                     fit_hill(assay = "endp", chem = "casn") |> 
+#'                     dplyr::filter(!tp.sd.imputed, !logAC50.sd.imputed)) |>
+#'   # Calculate response
+#'   calculate_response() |>
+#'   # Perform sensitivity analysis
+#'   sensitivity_analysis()
+#' 
+#' # Print GeoTox object
+#' geoTox
+#' 
+#' # Plot hill fits
+#' plot(geoTox, type = "hill")
+#' # Plot exposure data
+#' plot(geoTox, type = "exposure")
+#' # Plot response data
+#' plot(geoTox)
+#' plot(geoTox, assays = "TOX21_H2AX_HTRF_CHO_Agonist_ratio")
+#' # Plot sensitivity data
+#' plot(geoTox, type = "sensitivity")
+#' plot(geoTox, type = "sensitivity", assay = "TOX21_H2AX_HTRF_CHO_Agonist_ratio")
 GeoTox <- function() {
   structure(
     list(
@@ -169,7 +209,5 @@ plot.GeoTox <- function(x,
                   ncol            = ncol)
   } else if (type == "sensitivity") {
     plot_sensitivity(x, ...)
-  } else {
-    stop("Invalid type.")
   }
 }

R/plot_sensitivity.R

@@ -10,6 +10,25 @@
 #'
 #' @return ggplot2 object.
 #' @export
+#' 
+#' @examples
+#' # Required GeoTox fields are generated by first running [calculate_response]
+#' # and [sensitivity_analysis] on a GeoTox object. This will create the fields
+#' # `resp` and `sensitivity`. For this example, dummy data will be used.
+#' make_data <- function(n = 5, metric = "GCA.Eff") {
+#'   list(stats::setNames(data.frame(1:n, runif(n)),
+#'                        c("sample", metric)))
+#' }
+#' 
+#' x <- GeoTox()
+#' x$resp <- make_data()
+#' x$sensitivity <- list(age = make_data(),
+#'                       obesity = make_data(),
+#'                       css_params = make_data(),
+#'                       fit_params = make_data(),
+#'                       C_ext = make_data())
+#' 
+#' plot_sensitivity(x)
 plot_sensitivity <- function(x,
                              metric = "GCA.Eff",
                              assay = NULL,
@@ -20,6 +39,9 @@ plot_sensitivity <- function(x,
   if (is.null(x$sensitivity)) {
     stop("No sensitivity data found.", call. = FALSE)
   }
+  if (is.null(x$resp)) {
+    stop("No baseline response data found.", call. = FALSE)
+  }
 
   df <- get_sensitivity_df(x, metric = metric, assay = assay)
   fig <- plot_sensitivity_df(df, y = y, xlab = xlab, ylab = ylab)

_pkgdown.yml

@@ -18,6 +18,8 @@ reference:
   - compute_sensitivity
   - plot_sensitivity
   - sensitivity_analysis
+  - set_boundaries
+  - set_hill_params
   - simulate_population
 - subtitle: Simulate
   desc: Functions that perform Monte Carlo simulation
@@ -45,7 +47,6 @@ reference:
   - starts_with("get_")
   - starts_with("hill_")
   - starts_with("obj_")
-  - starts_with("set_")
 - subtitle: Datasets
   desc: Package datasets
 - contents:

tests/testthat/test-GeoTox.R

@@ -0,0 +1,112 @@
+test_that("errors", {
+
+  geoTox <- GeoTox()
+
+  expect_error(plot(geoTox),
+               "No response data found")
+
+  geoTox$resp <- list()
+
+  expect_error(plot(geoTox),
+               "No region boundary data found")
+})
+
+test_that("default params", {
+
+  geoTox <- GeoTox()
+
+  expect_equal(geoTox$par,
+               list(n = 1e3,
+                    IR_params = NULL,
+                    obesity = list(obes_prev  = "OBESITY_CrudePrev",
+                                   obes_sd    = "OBESITY_SD",
+                                   obes_label = "FIPS"),
+                    exposure = list(expos_mean  = "mean",
+                                    expos_sd    = "sd",
+                                    expos_label = "casn"),
+                    internal_dose = list(time    = 1,
+                                         BW      = 1,
+                                         scaling = 1),
+                    resp = list(max_mult = 1.5)))
+})
+
+test_that("package data subset", {
+
+  set.seed(2357)
+  n <- 10 # Population size
+  m <- 5 # Number of regions
+  idx <- if (m < 100) sample(1:100, m) else 1:100
+
+  expect_no_error(
+    geoTox <- GeoTox() |>
+      # Set region and group boundaries (for plotting)
+      set_boundaries(region = geo_tox_data$boundaries$county,
+                     group  = geo_tox_data$boundaries$state) |>
+      # Simulate populations for each region
+      simulate_population(age           = split(geo_tox_data$age, ~FIPS)[idx],
+                          obesity       = geo_tox_data$obesity[idx, ],
+                          exposure      = split(geo_tox_data$exposure, ~FIPS)[idx],
+                          simulated_css = geo_tox_data$simulated_css,
+                          n             = n) |> 
+      # Estimated Hill parameters
+      set_hill_params(geo_tox_data$dose_response |>
+                        fit_hill(assay = "endp", chem = "casn") |>
+                        dplyr::filter(!tp.sd.imputed, !logAC50.sd.imputed)) |>
+      # Calculate response
+      calculate_response() |>
+      # Perform sensitivity analysis
+      sensitivity_analysis()
+  )
+
+  # Print GeoTox object
+  expect_no_error(capture_output(print(geoTox)))
+
+  # Plot hill fits
+  expect_no_error(plot(geoTox, type = "hill"))
+  # Plot exposure data
+  expect_no_error(plot(geoTox, type = "exposure"))
+  # Plot response data
+  # Warning from resp_quantiles()
+  expect_warning(plot(geoTox),
+                 "Multiple assays found, using first assay")
+  expect_no_error(plot(geoTox, assays = "TOX21_H2AX_HTRF_CHO_Agonist_ratio"))
+  # Plot sensitivity data
+  expect_warning(expect_warning(
+    plot(geoTox, type = "sensitivity"),
+    "Multiple assays found, using first assay"),
+    "Removed \\d+ NA value")
+  expect_no_error(plot(geoTox,
+                       type = "sensitivity",
+                       assay = "TOX21_H2AX_HTRF_CHO_Agonist_ratio"))
+})
+
+test_that("print corner cases", {
+
+  geoTox <- GeoTox()
+
+  expect_no_error(capture_output(print(geoTox)))
+
+  geoTox$age <- list()
+  geoTox$obesity <- c(0)
+  geoTox$par <- NULL
+
+  expect_no_error(capture_output(print(geoTox)))
+
+  geoTox <- GeoTox() |>
+    set_hill_params(geo_tox_data$dose_response |>
+                      dplyr::filter(endp == "TOX21_H2AX_HTRF_CHO_Agonist_ratio",
+                                    casn == "510-15-6") |>
+                      fit_hill())
+
+  expect_no_error(capture_output(print(geoTox)))
+
+  geoTox <- geoTox |>
+    set_hill_params(geo_tox_data$dose_response |>
+                      fit_hill(assay = "endp", chem = "casn"))
+
+  expect_no_error(capture_output(print(geoTox)))
+
+  geoTox$resp <- c(1:5)
+
+  expect_no_error(capture_output(print(geoTox)))
+})

tests/testthat/test-plot_sensitivity.R

@@ -0,0 +1,89 @@
+test_that("errors", {
+
+  geoTox <- GeoTox()
+
+  expect_error(plot_sensitivity(geoTox),
+               "No sensitivity data found")
+
+  geoTox$sensitivity <- list()
+
+  expect_error(plot_sensitivity(geoTox),
+               "No baseline response data found")
+})
+
+test_that("no assay", {
+
+  make_data <- function(n = 5, metric = "GCA.Eff") {
+    list(stats::setNames(data.frame(1:n, runif(n)),
+                         c("sample", metric)))
+  }
+
+  x <- GeoTox()
+  n <- 5
+  x$resp <- make_data(n)
+  x$sensitivity <- list(age = make_data(n),
+                        obesity = make_data(n),
+                        css_params = make_data(n),
+                        fit_params = make_data(n),
+                        C_ext = make_data(n))
+
+  expect_no_error(plot_sensitivity(x))
+})
+
+test_that("with assay", {
+
+  make_data <- function(n = 5, metric = "GCA.Eff", assay = c("a1", "a2")) {
+    list(stats::setNames(data.frame(rep(assay, each = n),
+                                    rep(1:n, 2),
+                                    runif(2 * n)),
+                         c("assay", "sample", metric)))
+  }
+
+  x <- GeoTox()
+  n <- 5
+  x$resp <- make_data(n)
+  x$sensitivity <- list(age = make_data(n),
+                        obesity = make_data(n),
+                        css_params = make_data(n),
+                        fit_params = make_data(n),
+                        C_ext = make_data(n))
+
+  expect_warning(plot_sensitivity(x),
+                 "Multiple assays found, using first assay")
+  expect_no_error(plot_sensitivity(x, assay = "a1"))
+  expect_no_error(plot_sensitivity(x, assay = "a2"))
+})
+
+test_that("data NAs and missing", {
+
+  make_data <- function(n = 5, metric = "GCA.Eff", p_NA = 0.5) {
+    df <- stats::setNames(data.frame(1:n, runif(n)),
+                          c("sample", metric))
+    df[sample(n, n * p_NA), metric] <- NA
+    list(df)
+  }
+
+  x <- GeoTox()
+  n <- 5
+  x$resp <- make_data(n)
+  x$sensitivity <- list(age = make_data(n),
+                        obesity = make_data(n),
+                        css_params = make_data(n),
+                        fit_params = make_data(n),
+                        C_ext = make_data(n))
+
+  expect_warning(plot_sensitivity(x),
+                 "Removed \\d+ NA value")
+
+  p_NA <- 1
+  x$resp <- make_data(n, p_NA = p_NA)
+  x$sensitivity <- list(age = make_data(n, p_NA = p_NA),
+                        obesity = make_data(n, p_NA = p_NA),
+                        css_params = make_data(n, p_NA = p_NA),
+                        fit_params = make_data(n, p_NA = p_NA),
+                        C_ext = make_data(n, p_NA = p_NA))
+
+  expect_error(expect_warning(plot_sensitivity(x),
+                              "Removed \\d+ NA value"),
+               "No data to plot")
+})