diff --git a/DESCRIPTION b/DESCRIPTION
index f7d12c9..bcab1f4 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -30,7 +30,7 @@ Suggests:
knitr,
lubridate,
rmarkdown,
- testthat
+ testthat,
vdiffr
VignetteBuilder: knitr
LinkingTo:
diff --git a/NAMESPACE b/NAMESPACE
index 4a255a5..7fa16cf 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -1,5 +1,6 @@
# Generated by roxygen2: do not edit by hand
+S3method(plot,biokinetics_population_trajectories)
export(add_exposure_data)
export(biokinetics)
export(biokinetics_priors)
@@ -24,4 +25,5 @@ importFrom(ggplot2,geom_smooth)
importFrom(ggplot2,ggplot)
importFrom(ggplot2,guide_legend)
importFrom(ggplot2,guides)
+importFrom(ggplot2,scale_y_continuous)
useDynLib(epikinetics, .registration = TRUE)
diff --git a/R/biokinetics.R b/R/biokinetics.R
index eb3f7f3..3f6d48a 100644
--- a/R/biokinetics.R
+++ b/R/biokinetics.R
@@ -142,6 +142,11 @@ biokinetics <- R6::R6Class(
dt_out[, t_id := NULL]
+ if (private$scale == "natural") {
+ dt_out <- convert_log2_scale_inverse(
+ dt_out, vars_to_transform = "mu")
+ }
+
if (summarise == TRUE) {
logger::log_info("Summarising into quantiles")
dt_out <- summarise_draws(
@@ -266,6 +271,8 @@ biokinetics <- R6::R6Class(
)
},
#' @description Plot the kinetics trajectory predicted by the model priors.
+ #' Note that this is on a log scale, regardless of whether the data was provided
+ #' on a log or a natural scale.
#' @return A ggplot2 object.
#' @param tmax Integer. The number of time points in each simulated trajectory. Default 150.
#' @param n_draws Integer. The number of trajectories to simulate. Default 2000.
@@ -276,10 +283,12 @@ biokinetics <- R6::R6Class(
n_draws = n_draws,
data = private$data)
},
- #' @description Plot model input data with a smoothing function.
+ #' @description Plot model input data with a smoothing function. Note that
+ #' this plot is on a log scale, regardless of whether data was provided on a
+ #' log or a natural scale.
#' @return A ggplot2 object.
- plot_data = function() {
- plot_data(private$data)
+ plot_model_inputs = function() {
+ plot_data(private$data, private$all_formula_vars)
},
#' @description View the data that is passed to the stan model, for debugging purposes.
#' @return A list of arguments that will be passed to the stan model.
@@ -388,6 +397,7 @@ biokinetics <- R6::R6Class(
t_max = 150,
summarise = TRUE,
n_draws = 2500) {
+
private$check_fitted()
validate_numeric(t_max)
validate_logical(summarise)
@@ -404,17 +414,9 @@ biokinetics <- R6::R6Class(
dt_out <- dt_out[
, lapply(.SD, function(x) if (is.factor(x)) forcats::fct_drop(x) else x)]
- if (private$scale == "log") {
- return(dt_out)
- }
-
- if (summarise) {
- dt_out <- convert_log2_scale_inverse(
- dt_out, vars_to_transform = c("me", "lo", "hi"))
- } else {
- dt_out <- convert_log2_scale_inverse(
- dt_out, vars_to_transform = "mu")
- }
+ class(dt_out) <- append("biokinetics_population_trajectories", class(dt_out))
+ attr(dt_out, "summarised") <- summarise
+ attr(dt_out, "covariates") <- private$all_formula_vars
dt_out
},
#' @description Process the stan model results into a data.table.
@@ -550,7 +552,10 @@ biokinetics <- R6::R6Class(
by = c("calendar_day", "titre_type"))
}
- dt_out[, time_shift := time_shift]
+ dt_out <- dt_out[, time_shift := time_shift]
+ class(dt_out) <- append("biokinetics_individual_trajectories", class(dt_out))
+ attr(dt_out, "summarised") <- summarise
+ dt_out
}
)
)
diff --git a/R/epikinetics-package.R b/R/epikinetics-package.R
index 4404fc3..35b5299 100644
--- a/R/epikinetics-package.R
+++ b/R/epikinetics-package.R
@@ -11,7 +11,7 @@
#' @importFrom data.table .NGRP
#' @importFrom data.table .SD
#' @importFrom data.table data.table
-#' @importFrom ggplot2 aes facet_wrap geom_point geom_ribbon geom_line geom_smooth ggplot guides guide_legend
+#' @importFrom ggplot2 aes facet_wrap geom_point geom_ribbon geom_line geom_smooth ggplot guides guide_legend scale_y_continuous
#' @useDynLib epikinetics, .registration = TRUE
## usethis namespace: end
diff --git a/R/plot.R b/R/plot.R
index 3ebc65f..ed11e5b 100644
--- a/R/plot.R
+++ b/R/plot.R
@@ -7,22 +7,28 @@
#' @param priors A named list of type 'biokinetics_priors'.
#' @param tmax Integer. The number of time points in each simulated trajectory. Default 150.
#' @param n_draws Integer. The number of trajectories to simulate. Default 2000.
-#' @param data Optional data.frame with columns t_since_last_exp and value. The raw data to compare to.
+#' @param data Optional data.frame with columns time_since_last_exp and value. The raw data to compare to.
plot_prior_predictive <- function(priors,
tmax = 150,
n_draws = 2000,
data = NULL) {
validate_priors(priors)
+
+ # Declare variables to suppress notes when compiling package
+ # https://github.com/Rdatatable/data.table/issues/850#issuecomment-259466153
+ t0 <- tp <- ts <- m1 <- m2 <- m3 <- NULL
+ time_since_last_exp <- me <- lo <- hi <- value <- mu <- NULL
+
if (!is.null(data)) {
- validate_required_cols(data, c("t_since_last_exp", "value"))
+ validate_required_cols(data, c("time_since_last_exp", "value"))
}
params <- data.table(
- t0 = rnorm(n_draws, priors$mu_t0, priors$sigma_t0), # titre value at t0
- tp = rnorm(n_draws, priors$mu_tp, priors$sigma_tp), # time of peak
- ts = rnorm(n_draws, priors$mu_ts, priors$sigma_ts), # time of set point
- m1 = rnorm(n_draws, priors$mu_m1, priors$sigma_m1), # gradient 1
- m2 = rnorm(n_draws, priors$mu_m2, priors$sigma_m2), # gradient 2
- m3 = rnorm(n_draws, priors$mu_m3, priors$sigma_m3) # gradient 3
+ t0 = stats::rnorm(n_draws, priors$mu_t0, priors$sigma_t0), # titre value at t0
+ tp = stats::rnorm(n_draws, priors$mu_tp, priors$sigma_tp), # time of peak
+ ts = stats::rnorm(n_draws, priors$mu_ts, priors$sigma_ts), # time of set point
+ m1 = stats::rnorm(n_draws, priors$mu_m1, priors$sigma_m1), # gradient 1
+ m2 = stats::rnorm(n_draws, priors$mu_m2, priors$sigma_m2), # gradient 2
+ m3 = stats::rnorm(n_draws, priors$mu_m3, priors$sigma_m3) # gradient 3
)
times <- data.table(t = 1:tmax)
@@ -40,16 +46,66 @@ plot_prior_predictive <- function(priors,
geom_ribbon(aes(x = t, ymin = lo, ymax = hi), alpha = 0.5)
if (!is.null(data)) {
- plot <- plot + geom_point(data = data, aes(x = t_since_last_exp, y = value))
+ plot <- plot + geom_point(data = data, aes(x = time_since_last_exp, y = value))
}
plot
}
-plot_data <- function(data) {
- validate_required_cols(data, c("t_since_last_exp", "value", "titre_type"))
+plot_data <- function(data, covariates) {
+ validate_required_cols(data, c("time_since_last_exp", "value", "titre_type"))
+ # Declare variables to suppress notes when compiling package
+ # https://github.com/Rdatatable/data.table/issues/850#issuecomment-259466153
+ time_since_last_exp <- value <- titre_type <- NULL
+
ggplot(data) +
- geom_point(aes(x = t_since_last_exp, y = value, colour = titre_type)) +
- geom_smooth(aes(x = t_since_last_exp, y = value, colour = titre_type)) +
- facet_wrap(~titre_type) +
+ geom_point(aes(x = time_since_last_exp, y = value, colour = titre_type)) +
+ geom_smooth(aes(x = time_since_last_exp, y = value, colour = titre_type)) +
+ facet_wrap(eval(parse(text = facet_formula(covariates)))) +
guides(colour = guide_legend(title = "Titre type"))
}
+
+#' Plot method for "biokinetics_population_trajectories" class
+#'
+#' @param x An object of class "biokinetics_population_trajectories". These are
+#' generated by running biokinetics$simulate_populate_trajectories(). See
+#' \href{../../epikinetics/html/biokinetics.html#method-biokinetics-simulate_population_trajectories}{\code{biokinetics$simulate_population_trajectories()}}
+#' @param \dots Further arguments passed to the method. For now the only extra argument supported is "data", which should be
+#' a data.table containing raw data as provided to the biokinetics model.
+#' @export
+plot.biokinetics_population_trajectories <- function(x, ...) {
+ if (!attr(x, "summarised")) {
+ by <- setdiff(colnames(x), c("t0_pop", "tp_pop", "ts_pop",
+ "m1_pop", "m2_pop", "m3_pop",
+ ".draw", "mu"))
+ x <- summarise_draws(
+ x, column_name = "mu", by = by)
+ }
+ args <- list(...)
+ # Declare variables to suppress notes when compiling package
+ # https://github.com/Rdatatable/data.table/issues/850#issuecomment-259466153
+ time_since_last_exp <- value <- me <- titre_type <- lo <- hi <- NULL
+ day <- last_exp_day <- NULL
+
+ plot <- ggplot(x) +
+ geom_line(aes(x = time_since_last_exp, y = me, colour = titre_type)) +
+ geom_ribbon(aes(x = time_since_last_exp,
+ ymin = lo,
+ ymax = hi,
+ fill = titre_type), alpha = 0.5)
+
+ if (!is.null(args$data)) {
+ data <- args$data
+ validate_required_cols(data)
+ plot <- plot +
+ geom_point(data = data,
+ aes(x = as.integer(day - last_exp_day, units = "days"),
+ y = value), size = 0.5, alpha = 0.5)
+ }
+ plot + scale_y_continuous(trans = "log2") +
+ facet_wrap(eval(parse(text = facet_formula(attr(x, "covariates"))))) +
+ guides(fill = guide_legend(title = "Titre type"), colour = "none")
+}
+
+facet_formula <- function(covariates) {
+ paste("~", paste(c("titre_type", covariates), collapse = "+"))
+}
diff --git a/R/utils.R b/R/utils.R
index 1ee1a2e..6d770f4 100644
--- a/R/utils.R
+++ b/R/utils.R
@@ -1,5 +1,5 @@
convert_log2_scale <- function(
- dt_in, vars_to_transform = "titre",
+ dt_in, vars_to_transform = "value",
simplify_limits = TRUE) {
dt_out <- data.table::copy(dt_in)
diff --git a/man/biokinetics.Rd b/man/biokinetics.Rd
index 842b49f..feb9501 100644
--- a/man/biokinetics.Rd
+++ b/man/biokinetics.Rd
@@ -12,7 +12,7 @@ fit it to a dataset.
\itemize{
\item \href{#method-biokinetics-new}{\code{biokinetics$new()}}
\item \href{#method-biokinetics-plot_prior_predictive}{\code{biokinetics$plot_prior_predictive()}}
-\item \href{#method-biokinetics-plot_data}{\code{biokinetics$plot_data()}}
+\item \href{#method-biokinetics-plot_model_inputs}{\code{biokinetics$plot_model_inputs()}}
\item \href{#method-biokinetics-get_stan_data}{\code{biokinetics$get_stan_data()}}
\item \href{#method-biokinetics-get_covariate_lookup_table}{\code{biokinetics$get_covariate_lookup_table()}}
\item \href{#method-biokinetics-fit}{\code{biokinetics$fit()}}
@@ -68,6 +68,8 @@ An epikinetics::biokinetics object.
\if{latex}{\out{\hypertarget{method-biokinetics-plot_prior_predictive}{}}}
\subsection{Method \code{plot_prior_predictive()}}{
Plot the kinetics trajectory predicted by the model priors.
+Note that this is on a log scale, regardless of whether the data was provided
+on a log or a natural scale.
\subsection{Usage}{
\if{html}{\out{
}}\preformatted{biokinetics$plot_prior_predictive(tmax = 150, n_draws = 2000)}\if{html}{\out{
}}
}
@@ -86,12 +88,14 @@ A ggplot2 object.
}
}
\if{html}{\out{ }}
-\if{html}{\out{ }}
-\if{latex}{\out{\hypertarget{method-biokinetics-plot_data}{}}}
-\subsection{Method \code{plot_data()}}{
-Plot model input data with a smoothing function.
+\if{html}{\out{ }}
+\if{latex}{\out{\hypertarget{method-biokinetics-plot_model_inputs}{}}}
+\subsection{Method \code{plot_model_inputs()}}{
+Plot model input data with a smoothing function. Note that
+this plot is on a log scale, regardless of whether data was provided on a
+log or a natural scale.
\subsection{Usage}{
-\if{html}{\out{}}\preformatted{biokinetics$plot_data()}\if{html}{\out{
}}
+\if{html}{\out{}}\preformatted{biokinetics$plot_model_inputs()}\if{html}{\out{
}}
}
\subsection{Returns}{
diff --git a/man/plot.biokinetics_population_trajectories.Rd b/man/plot.biokinetics_population_trajectories.Rd
new file mode 100644
index 0000000..aeecc78
--- /dev/null
+++ b/man/plot.biokinetics_population_trajectories.Rd
@@ -0,0 +1,19 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/plot.R
+\name{plot.biokinetics_population_trajectories}
+\alias{plot.biokinetics_population_trajectories}
+\title{Plot method for "biokinetics_population_trajectories" class}
+\usage{
+\method{plot}{biokinetics_population_trajectories}(x, ...)
+}
+\arguments{
+\item{x}{An object of class "biokinetics_population_trajectories". These are
+generated by running biokinetics$simulate_populate_trajectories(). See
+\href{../../epikinetics/html/biokinetics.html#method-biokinetics-simulate_population_trajectories}{\code{biokinetics$simulate_population_trajectories()}}}
+
+\item{\dots}{Further arguments passed to the method. For now the only extra argument supported is "data", which should be
+a data.table containing raw data as provided to the biokinetics model.}
+}
+\description{
+Plot method for "biokinetics_population_trajectories" class
+}
diff --git a/man/plot_prior_predictive.Rd b/man/plot_prior_predictive.Rd
index c1f64f7..fc6aa0c 100644
--- a/man/plot_prior_predictive.Rd
+++ b/man/plot_prior_predictive.Rd
@@ -14,7 +14,7 @@ plot_prior_predictive(priors, tmax = 150, n_draws = 2000, data = NULL)
\item{n_draws}{Integer. The number of trajectories to simulate. Default 2000.}
-\item{data}{Optional data.frame with columns t_since_last_exp and value. The raw data to compare to.}
+\item{data}{Optional data.frame with columns time_since_last_exp and value. The raw data to compare to.}
}
\value{
A ggplot2 object.
diff --git a/tests/testthat/_snaps/plots/inputdata-covariates.svg b/tests/testthat/_snaps/plots/inputdata-covariates.svg
new file mode 100644
index 0000000..20b5d78
--- /dev/null
+++ b/tests/testthat/_snaps/plots/inputdata-covariates.svg
@@ -0,0 +1,2506 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+Ancestral
+
+Previously infected (Pre-Omicron)
+
+
+
+
+
+
+
+
+
+
+Delta
+
+Infection naive
+
+
+
+
+
+
+
+
+
+
+Delta
+
+Previously infected (Pre-Omicron)
+
+
+
+
+
+
+
+
+
+
+Alpha
+
+Infection naive
+
+
+
+
+
+
+
+
+
+
+Alpha
+
+Previously infected (Pre-Omicron)
+
+
+
+
+
+
+
+
+
+
+Ancestral
+
+Infection naive
+
+
+
+
+
+
+0
+200
+400
+600
+
+
+
+
+0
+200
+400
+600
+
+
+
+
+0
+200
+400
+600
+0
+3
+6
+9
+
+
+
+
+0
+3
+6
+9
+
+
+
+
+time_since_last_exp
+value
+
+Titre type
+
+
+
+
+
+
+
+
+
+
+
+
+Alpha
+Ancestral
+Delta
+inputdata_covariates
+
+
diff --git a/tests/testthat/_snaps/plots/inputdata.svg b/tests/testthat/_snaps/plots/inputdata.svg
index c0b05fd..e4c509f 100644
--- a/tests/testthat/_snaps/plots/inputdata.svg
+++ b/tests/testthat/_snaps/plots/inputdata.svg
@@ -2384,7 +2384,7 @@
-t_since_last_exp
+time_since_last_exp
value
Titre type
diff --git a/tests/testthat/_snaps/plots/populationtrajectories-data.svg b/tests/testthat/_snaps/plots/populationtrajectories-data.svg
new file mode 100644
index 0000000..6340bda
--- /dev/null
+++ b/tests/testthat/_snaps/plots/populationtrajectories-data.svg
@@ -0,0 +1,2500 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+Ancestral
+
+Previously infected (Pre-Omicron)
+
+
+
+
+
+
+
+
+
+
+Delta
+
+Infection naive
+
+
+
+
+
+
+
+
+
+
+Delta
+
+Previously infected (Pre-Omicron)
+
+
+
+
+
+
+
+
+
+
+Alpha
+
+Infection naive
+
+
+
+
+
+
+
+
+
+
+Alpha
+
+Previously infected (Pre-Omicron)
+
+
+
+
+
+
+
+
+
+
+Ancestral
+
+Infection naive
+
+
+
+
+
+
+0
+200
+400
+600
+
+
+
+
+0
+200
+400
+600
+
+
+
+
+0
+200
+400
+600
+128
+8192
+524288
+33554432
+
+
+
+
+128
+8192
+524288
+33554432
+
+
+
+
+time_since_last_exp
+me
+
+Titre type
+
+
+
+
+
+
+Alpha
+Ancestral
+Delta
+populationtrajectories_data
+
+
diff --git a/tests/testthat/_snaps/plots/populationtrajectories.svg b/tests/testthat/_snaps/plots/populationtrajectories.svg
new file mode 100644
index 0000000..3c2d3b9
--- /dev/null
+++ b/tests/testthat/_snaps/plots/populationtrajectories.svg
@@ -0,0 +1,245 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+Ancestral
+
+Previously infected (Pre-Omicron)
+
+
+
+
+
+
+
+
+
+
+Delta
+
+Infection naive
+
+
+
+
+
+
+
+
+
+
+Delta
+
+Previously infected (Pre-Omicron)
+
+
+
+
+
+
+
+
+
+
+Alpha
+
+Infection naive
+
+
+
+
+
+
+
+
+
+
+Alpha
+
+Previously infected (Pre-Omicron)
+
+
+
+
+
+
+
+
+
+
+Ancestral
+
+Infection naive
+
+
+
+
+
+
+0
+50
+100
+150
+
+
+
+
+0
+50
+100
+150
+
+
+
+
+0
+50
+100
+150
+128
+8192
+524288
+33554432
+
+
+
+
+128
+8192
+524288
+33554432
+
+
+
+
+time_since_last_exp
+me
+
+Titre type
+
+
+
+
+
+
+Alpha
+Ancestral
+Delta
+populationtrajectories
+
+
diff --git a/tests/testthat/test-data.R b/tests/testthat/test-data.R
index 57ebbaf..f1ccbd1 100644
--- a/tests/testthat/test-data.R
+++ b/tests/testthat/test-data.R
@@ -52,7 +52,9 @@ test_that("Natural scale data is converted to log scale for stan", {
dat <- data.table::fread(system.file("delta_full.rds", package = "epikinetics"))
mod <- biokinetics$new(data = dat)
stan_data <- mod$get_stan_data()
- expect_equal(stan_data$value, convert_log2_scale(dat, "value")$value, ignore_attr = TRUE)
+ expect_equal(stan_data$value,
+ convert_log2_scale(dat, "value")$value,
+ ignore_attr = TRUE)
})
test_that("Log scale data is passed directly to stan", {
diff --git a/tests/testthat/test-plots.R b/tests/testthat/test-plots.R
index d01a8e1..32723b3 100644
--- a/tests/testthat/test-plots.R
+++ b/tests/testthat/test-plots.R
@@ -8,8 +8,8 @@ test_that("Can plot prior prediction up to tmax", {
test_that("Can plot prior prediction with data points", {
data <- data.table::fread(system.file("delta_full.rds", package = "epikinetics"))
priors <- biokinetics_priors()
- expect_error(plot_prior_predictive(priors, data = data), "Missing required columns: t_since_last_exp")
- data[, `:=`(t_since_last_exp = as.integer(day - last_exp_day, units = "days"))]
+ expect_error(plot_prior_predictive(priors, data = data), "Missing required columns: time_since_last_exp")
+ data[, `:=`(time_since_last_exp = as.integer(day - last_exp_day, units = "days"))]
plot <- plot_prior_predictive(priors, data = data, n_draws = 500)
expect_equal(length(plot$layers), 3)
})
@@ -28,7 +28,6 @@ test_that("Can plot prior predictions from model", {
})
test_that("Prior predictions from model are the same", {
- skip_on_ci()
data <- data.table::fread(system.file("delta_full.rds", package = "epikinetics"))
priors <- biokinetics_priors(mu_values = c(4.1, 11, 65, 0.2, -0.01, 0.01),
sigma_values = c(2.0, 2.0, 3.0, 0.01, 0.01, 0.001))
@@ -40,10 +39,45 @@ test_that("Prior predictions from model are the same", {
vdiffr::expect_doppelganger("priorpredictive", plot)
})
-test_that("Plotted data is are the same", {
- skip_on_ci()
+test_that("Can plot input data", {
data <- data.table::fread(system.file("delta_full.rds", package = "epikinetics"))
mod <- biokinetics$new(data = data)
- plot <- mod$plot_data()
+ plot <- mod$plot_model_inputs()
vdiffr::expect_doppelganger("inputdata", plot)
+
+ mod <- biokinetics$new(data = data, covariate_formula = ~0 + infection_history)
+ plot <- mod$plot_model_inputs()
+ vdiffr::expect_doppelganger("inputdata_covariates", plot)
+})
+
+mock_model <- function(name, package) {
+ list(sample = function(x, ...) readRDS(test_path("testdata", "testdraws.rds")))
+}
+
+test_that("Summarised and un-summarised population trajectories give same plots", {
+ # note that this is using a pre-fitted model with very few iterations, so the
+ # fits won't look very good
+ local_mocked_bindings(
+ stan_package_model = mock_model, .package = "instantiate"
+ )
+ mod <- biokinetics$new(file_path = system.file("delta_full.rds", package = "epikinetics"),
+ covariate_formula = ~0 + infection_history)
+ mod$fit()
+ trajectories <- mod$simulate_population_trajectories(summarise = TRUE)
+ unsummarised_trajectories <- mod$simulate_population_trajectories(summarise = FALSE)
+ vdiffr::expect_doppelganger("populationtrajectories", plot(trajectories))
+ vdiffr::expect_doppelganger("populationtrajectories", plot(unsummarised_trajectories))
+})
+
+test_that("Can plot population trajectories with data", {
+ local_mocked_bindings(
+ stan_package_model = mock_model, .package = "instantiate"
+ )
+ # note that this is using a pre-fitted model with very few iterations, so the
+ # fits won't look very good
+ data <- data.table::fread(system.file("delta_full.rds", package = "epikinetics"))
+ mod <- biokinetics$new(data = data, covariate_formula = ~0 + infection_history)
+ fit <- mod$fit()
+ trajectories <- mod$simulate_population_trajectories()
+ vdiffr::expect_doppelganger("populationtrajectories_data", plot(trajectories, data = data))
})