diff --git a/R/adherence.R b/R/new_adherence.R similarity index 100% rename from R/adherence.R rename to R/new_adherence.R diff --git a/tests/testthat/setup.R b/tests/testthat/setup.R index 906863f0..1e9878ac 100644 --- a/tests/testthat/setup.R +++ b/tests/testthat/setup.R @@ -1,6 +1,13 @@ +# ============================================================================= +# Shared models compiled once at test session startup +# ============================================================================= + +# --- Library models (standard PK models) --- mod_1cmt_iv <- new_ode_model("pk_1cmt_iv") mod_2cmt_iv <- new_ode_model("pk_2cmt_iv") mod_1cmt_oral <- new_ode_model("pk_1cmt_oral") + +# --- Custom 1cmt oral with lagtime --- mod_1cmt_oral_lagtime <- new_ode_model( code = " dAdt[0] = -KA * A[0] @@ -11,7 +18,9 @@ mod_1cmt_oral_lagtime <- new_ode_model( dose = list(cmt = 1, bioav = 1), parameters = list(CL = 5, V = 50, KA = 0.5, TLAG = 0.83) ) -oral_1cmt_allometric <- new_ode_model( # also timevarying and dose-dependence factor + +# --- 1cmt oral with allometric scaling, time-varying CL, dose-dependent V --- +oral_1cmt_allometric <- new_ode_model( code = " if(t<168.0) { CLi = CL * pow(WT/70, 0.75) @@ -32,3 +41,167 @@ oral_1cmt_allometric <- new_ode_model( # also timevarying and dose-dependence fa declare_variables = c("CLi", "Vi"), parameters = list(KA = 0.5, CL = 5, V = 50) ) + +# --- From test_iov.R: IOV models --- +pars_iov <- list( + "kappa_CL_1" = 0, + "kappa_CL_2" = 0, + "kappa_CL_3" = 0, + "eta_CL" = 0, + "CL" = 5, + "V" = 50, + "KA" = 1 +) +pars_iov_no_iov <- list( + "CL" = 5, + "V" = 50, + "KA" = 1 +) +pk_iov_none <- new_ode_model( + code = " + dAdt[1] = -KA * A[1] + dAdt[2] = +KA * A[1] -(CL/V) * A[2] + ", + obs = list(cmt = 2, scale = "V"), + dose = list(cmt = 1, bioav = 1), + parameters = names(pars_iov_no_iov), + cpp_show_code = FALSE +) +pk_iov <- new_ode_model( + code = " + CL_iov = CL * exp(kappa_CL + eta_CL); + dAdt[1] = -KA * A[1] + dAdt[2] = +KA * A[1] -(CL_iov/V) * A[2] + ", + iov = list( + cv = list(CL = 0.2), + n_bins = 3 + ), + obs = list(cmt = 2, scale = "V"), + dose = list(cmt = 1, bioav = 1), + declare_variables = c("kappa_CL", "CL_iov"), + parameters = names(pars_iov), + cpp_show_code = FALSE +) + +# --- From test_compare_results.R --- +dose_in_cmt_2 <- new_ode_model( + code = " + dAdt[1] = -KA * A[1]; + dAdt[2] = KA*A[1] -(CL/V) * A[2] + dAdt[3] = S2*(A[2]-A[3]) + ", + obs = list(cmt = 2, scale = "V"), + dose = list(cmt = 2), + cpp_show_code = FALSE +) + +# --- From test_multi_obs.R: Multi-observation model --- +vars_multi_obs <- c("CONC", "METAB", "METAB2", "ACT") +pk_multi_obs <- new_ode_model( + code = "dAdt[1] = -(CL/V)*A[1]; CONC = 1000*A[1]/V; METAB = CONC/2; METAB2 = CONC * t; ACT = 15", + obs = list(variable = vars_multi_obs), + declare_variables = vars_multi_obs, + cpp_show_code = FALSE +) + +# --- From test_mixture_model.R --- +covs_mixture <- list(WT = new_covariate(70)) +mod_mixture <- new_ode_model( + code = " + dAdt[0] = -(CL*(WT/70.0)/V)*A[0]; + ", + pk_code = " ", + obs = list(cmt = 1, scale = "V"), + mixture = list(CL = list(values = c(5, 15), probability = 0.3)), + covariates = covs_mixture +) + +# --- Conditional models (skip on CRAN due to compilation time) --- +if (identical(Sys.getenv("NOT_CRAN"), "true")) { + + # Library models + pk_3cmt_iv <- new_ode_model("pk_3cmt_iv") + pk_2cmt_oral <- new_ode_model("pk_2cmt_oral") + pk_3cmt_oral <- new_ode_model("pk_3cmt_oral") + mod_1cmt_iv_mm <- new_ode_model("pk_1cmt_iv_mm") + + # --- From test_advan_with_auc.R: Models with AUC compartments --- + parameters_advan_auc <- list( + CL = 10, V = 50, KA = 0.5, Q = 5, V2 = 100, Q2 = 3, V3 = 150, F1 = 1 + ) + mod_1cmt_auc <- new_ode_model( + code = "dAdt[1] = -(CL/V)*A[1]; dAdt[2] = A[1]/V;", + parameters = parameters_advan_auc + ) + mod_2cmt_auc <- new_ode_model( + code = " + dAdt[1] = -(CL/V)*A[1] - (Q/V)*A[1] + (Q/V2)*A[2]; + dAdt[2] = +(Q/V)*A[1] - (Q/V2)*A[2]; + dAdt[3] = A[1]/V; + ", + parameters = parameters_advan_auc + ) + mod_3cmt_auc <- new_ode_model( + code = " + dAdt[1] = -(CL/V)*A[1] - (Q/V)*A[1] + (Q/V2)*A[2] - (Q2/V)*A[1] + (Q2/V3)*A[3]; + dAdt[2] = (Q/V)*A[1] -(Q/V2)*A[2] ; + dAdt[3] = (Q2/V)*A[1] - (Q2/V3)*A[3]; + dAdt[4] = A[1]/V; + ", + parameters = parameters_advan_auc + ) + + # --- From test_timevar_cov.R: 2cmt with time-varying covariates --- + mod_2cmt_timevar <- new_ode_model( + code = " + dAdt[1] = -(Q/V)*A[1] + (Q/V2)*A[2] -(CLi/V)*A[1]; + dAdt[2] = -(Q/V2)*A[2] + (Q/V)*A[1]; + ", + pk_code = "CLi = CL + CRCL", + obs = list(cmt = 2, scale = "V"), + covariates = list(CRCL = new_covariate(5)), + declare_variables = "CLi", + cpp = FALSE + ) + + # --- From test_t_init.R: Model with state initialization --- + mod_t_init <- new_ode_model( + code = "CLi = CL; Vi = V; dAdt[1] = -(CLi/Vi)*A[1]; CONC = A[1]/Vi", + state_init = "A[1] = TDM_INIT * Vi", + parameters = list(CL = 7.67, V = 97.7, TDM_INIT = 500), + obs = list(cmt = 1, scale = "Vi"), + declare_variables = c("CONC", "CLi", "Vi"), + cpp_show_code = FALSE + ) + + # --- From test_reparametrization.R: Carreno 2cmt model --- + covs_carreno <- list( + CRCL = new_covariate(5), + CL_HEMO = new_covariate(0) + ) + model_carreno <- new_ode_model( + code = " + CLi = SCLInter + SCLSlope * (CRCL*16.667) + CL_HEMO \ + Vi = V \ + Qi = K12 * Vi \ + V2i = Qi / K21 \ + dAdt[0] = -(CLi/V)*A[0] - K12*A[0] + K21*A[1] \ + dAdt[1] = K12*A[0] - K21*A[1] \ + dAdt[2] = A[0]/V + ", + parameters = list( + V = 25.76, SCLSlope = 0.036, K12 = 2.29, K21 = 1.44, + SCLInter = 0.18, TDM_INIT = 0 + ), + declare_variables = c("CLi", "Qi", "Vi", "V2i"), + covariates = covs_carreno, + obs = list(cmt = 1, scale = "V"), + reparametrization = list( + "CL" = "SCLInter + SCLSlope * (CRCL*16.667) + CL_HEMO", + "V" = "V", + "Q" = "K12 * V", + "V2" = "(K12 * V) / K21" + ) + ) +} diff --git a/tests/testthat/test_advan.R b/tests/testthat/test_advan.R index 5a8dce3a..c5bc8b25 100644 --- a/tests/testthat/test_advan.R +++ b/tests/testthat/test_advan.R @@ -124,3 +124,213 @@ test_that("Three compartment IV oral", { expect_true(!any(is.na(res3_oral$DV))) expect_equal(res3_oral, res3_oral_c) }) + +describe("ADVANs with AUC", { + # Uses models and parameters defined in setup.R (conditional, NOT_CRAN only): + # - mod_1cmt_auc, mod_2cmt_auc, mod_3cmt_auc + # - parameters_advan_auc + + if (identical(Sys.getenv("NOT_CRAN"), "true")) { + dose <- 100 + interval <- 12 + n_days <- 5 + t_inf <- 1.5 + t_obs <- c(3, 6, 8, 23, 47) + + ## bolus dataset + reg_bolus <- new_regimen( + amt = dose, + times = seq(0, interval * n_days * (24/interval), interval), + t_inf = t_inf, + type = "bolus" + ) + data_bolus <- advan_create_data( + reg_bolus, + parameters = parameters_advan_auc, + cmts = 5, + t_obs = t_obs + ) + + ## Infusion dataset + reg_infusion <- new_regimen( + amt = dose, + times = seq(0, interval * n_days * (24/interval), interval), + t_inf = t_inf, + type = "infusion" + ) + data_infusion <- advan_create_data( + reg_infusion, + parameters = parameters_advan_auc, + cmts = 6, + t_obs = t_obs + ) + } + + test_that("One compartment bolus ADVAN runs", { + skip_on_cran() + res1_iv_r <- advan("1cmt_iv_bolus", cpp=FALSE)(data_bolus) + res1_iv_c <- advan("1cmt_iv_bolus", cpp=TRUE)(data_bolus) + res1_iv_ode <- sim(ode = mod_1cmt_auc, regimen = reg_bolus, parameters = parameters_advan_auc, t_obs = t_obs) + + # AUC R + expect_equal(round(res1_iv_r[res1_iv_r$TIME %in% t_obs,]$AUC, 5), round(res1_iv_ode[res1_iv_ode$comp == 2,]$y, 5)) + + #AUC-C + expect_equal(round(res1_iv_c[res1_iv_c$TIME %in% t_obs,]$AUC, 5), round(res1_iv_ode[res1_iv_ode$comp == 2,]$y, 5)) + }) + + test_that("Two compartment bolus ADVAN runs", { + skip_on_cran() + res2_iv_r <- advan("2cmt_iv_bolus", cpp=FALSE)(data_bolus) + res2_iv_c <- advan("2cmt_iv_bolus", cpp=TRUE)(data_bolus) + res2_iv_ode <- sim(ode = mod_2cmt_auc, regimen = reg_bolus, parameters = parameters_advan_auc, t_obs = t_obs) + expect_equal( + round(res2_iv_r[res2_iv_r$TIME %in% t_obs,]$AUC, 5), + round(res2_iv_c[res2_iv_c$TIME %in% t_obs,]$AUC, 5) + ) + # AUC R + expect_equal( + round(res2_iv_r[res2_iv_r$TIME %in% t_obs,]$AUC, 5), + round(res2_iv_ode[res2_iv_ode$comp == 3,]$y, 5) + ) + + #AUC-C + expect_equal( + round(res2_iv_c[res2_iv_c$TIME %in% t_obs,]$AUC, 5), + round(res2_iv_ode[res2_iv_ode$comp == 3,]$y, 5) + ) + }) + + test_that("Two compartment infusion ADVAN runs", { + skip_on_cran() + res2_inf_r <- advan("2cmt_iv_infusion", cpp=FALSE)(data_infusion) + res2_inf_c <- advan("2cmt_iv_infusion", cpp=TRUE)(data_infusion) + res2_inf_ode <- sim(ode = mod_2cmt_auc, regimen = reg_infusion, parameters = parameters_advan_auc, t_obs = t_obs) + + expect_equal( + round(res2_inf_r[res2_inf_r$TIME %in% t_obs,]$AUC, 5), + round(res2_inf_c[res2_inf_c$TIME %in% t_obs,]$AUC, 5) + ) + + # AUC R + expect_equal( + round(res2_inf_r[res2_inf_r$TIME %in% t_obs,]$AUC, 5), + round(res2_inf_ode[res2_inf_ode$comp == 3,]$y, 5) + ) + + #AUC-C + expect_equal( + round(res2_inf_c[res2_inf_c$TIME %in% t_obs,]$AUC, 5), + round(res2_inf_ode[res2_inf_ode$comp == 3,]$y, 5) + ) + + }) + + test_that("Three compartment bolus ADVAN runs", { + skip_on_cran() + res3_iv_r <- advan("3cmt_iv_bolus", cpp=FALSE)(data_bolus) + res3_iv_c <- advan("3cmt_iv_bolus", cpp=TRUE)(data_bolus) + res3_iv_ode <- sim(ode = mod_3cmt_auc, regimen = reg_bolus, parameters = parameters_advan_auc, t_obs = t_obs) + expect_equal( + round(res3_iv_r[res3_iv_r$TIME %in% t_obs,]$AUC, 5), + round(res3_iv_c[res3_iv_c$TIME %in% t_obs,]$AUC, 5) + ) + # AUC R + expect_equal( + round(res3_iv_r[res3_iv_r$TIME %in% t_obs,]$AUC, 5), + round(res3_iv_ode[res3_iv_ode$comp == 4,]$y, 5) + ) + + #AUC-C + expect_equal( + round(res3_iv_c[res3_iv_c$TIME %in% t_obs,]$AUC, 5), + round(res3_iv_ode[res3_iv_ode$comp == 4,]$y, 5) + ) + }) + + test_that("Three compartment iv ADVAN runs", { + skip_on_cran() + res3_iv_r <- advan("3cmt_iv_infusion", cpp=FALSE)(data_infusion) + res3_iv_c <- advan("3cmt_iv_infusion", cpp=TRUE)(data_infusion) + res3_iv_ode <- sim(ode = mod_3cmt_auc, regimen = reg_infusion, parameters = parameters_advan_auc, t_obs = t_obs) + expect_equal( + round(res3_iv_r[res3_iv_r$TIME %in% t_obs,]$AUC, 5), + round(res3_iv_c[res3_iv_c$TIME %in% t_obs,]$AUC, 5) + ) + # AUC R + expect_equal( + round(res3_iv_r[res3_iv_r$TIME %in% t_obs,]$AUC, 5), + round(res3_iv_ode[res3_iv_ode$comp == 4,]$y, 5) + ) + + #AUC-C + expect_equal( + round(res3_iv_c[res3_iv_c$TIME %in% t_obs,]$AUC, 5), + round(res3_iv_ode[res3_iv_ode$comp == 4,]$y, 5) + ) + }) + +}) + +describe("ADVANs with covariates", { + test_that("Analytic and ODE models with covariates are the same", { + skip_on_cran() + + ## Create dataset + dose <- 100 + interval <- 12 + n_days <- 2 + parameters <- list( + CL = 10, + V = 50, + KA = 0.5, + Q = 5, + V2 = 100, + Q2 = 3, + V3 = 150, + F1 = 1 + ) + t_obs <- seq(0, 40, .1) + reg_bolus <- new_regimen( + amt = dose, + times = seq(0, interval * n_days * (24/interval), interval), + type = "bolus" + ) + ## TODO: there is slight difference in how bolus doses are handled. + ## Analytical equation is perhaps more consistent, so not testing + ## simulations at dose times. Should look into later. + t_obs <- t_obs[! t_obs %in% reg_bolus$dose_times] + covariates <- list(WT = new_covariate(80), CRCL=new_covariate(4.5)) + + ## Using analytic equations model: + data_ana <- sim( + analytical = "1cmt_iv_bolus", + parameters = parameters, + covariates = covariates, + regimen = reg_bolus, + t_obs = t_obs, + covariate_model = "CL = CL * (CRCL / 3)^0.75; V = V * (WT / 70.0)" + ) + + ## Using ODE model: + mod1 <- new_ode_model( + code = " + dAdt[1] = -( (CL*pow(CRCL/3.0, 0.75)) / (V*WT/70.0) ) * A[1]; + ", + covariates = covariates, + obs = list(cmt = 1, scale = "V*WT/70.0"), dose = list(cmt = 1) + ) + data_ode <- sim( + ode = mod1, + parameters = parameters, + covariates = covariates, + regimen = reg_bolus, + t_obs = t_obs, + duplicate_t_obs = TRUE, + only_obs = TRUE + ) + + expect_equal(nrow(data_ana), nrow(data_ode)) + expect_equal(round(data_ana$y,4), round(data_ode$y, 4)) + }) +}) \ No newline at end of file diff --git a/tests/testthat/test_advan_with_auc.R b/tests/testthat/test_advan_with_auc.R deleted file mode 100644 index ced7ef98..00000000 --- a/tests/testthat/test_advan_with_auc.R +++ /dev/null @@ -1,173 +0,0 @@ -if (identical(Sys.getenv("NOT_CRAN"), "true")) { - dose <- 100 - interval <- 12 - n_days <- 5 - t_inf <- 1.5 - parameters <- list( - CL = 10, - V = 50, - KA = 0.5, - Q = 5, - V2 = 100, - Q2 = 3, - V3 = 150, - F1 = 1 - ) - t_obs <- c(3, 6, 8, 23, 47) - - ## ODE models for testing - mod_1cmt <- new_ode_model( - code="dAdt[1] = -(CL/V)*A[1]; dAdt[2] = A[1]/V;", - parameters = parameters - ) - mod_2cmt <- new_ode_model( - code=" - dAdt[1] = -(CL/V)*A[1] - (Q/V)*A[1] + (Q/V2)*A[2]; - dAdt[2] = +(Q/V)*A[1] - (Q/V2)*A[2]; - dAdt[3] = A[1]/V; - ", - parameters = parameters - ) - mod_3cmt <- new_ode_model( - code=" - dAdt[1] = -(CL/V)*A[1] - (Q/V)*A[1] + (Q/V2)*A[2] - (Q2/V)*A[1] + (Q2/V3)*A[3]; - dAdt[2] = (Q/V)*A[1] -(Q/V2)*A[2] ; - dAdt[3] = (Q2/V)*A[1] - (Q2/V3)*A[3]; - dAdt[4] = A[1]/V; - ", - parameters = parameters - ) - - ## bolus dataset - reg_bolus <- new_regimen( - amt = dose, - times = seq(0, interval * n_days * (24/interval), interval), - t_inf = t_inf, - type = "bolus" - ) - data_bolus <- advan_create_data( - reg_bolus, - parameters = parameters, - cmts = 5, - t_obs = t_obs - ) - - ## Infusion dataset - reg_infusion <- new_regimen( - amt = dose, - times = seq(0, interval * n_days * (24/interval), interval), - t_inf = t_inf, - type = "infusion" - ) - data_infusion <- advan_create_data( - reg_infusion, - parameters = parameters, - cmts = 6, - t_obs = t_obs - ) -} - -test_that("One compartment bolus ADVAN runs", { - skip_on_cran() - res1_iv_r <- advan("1cmt_iv_bolus", cpp=FALSE)(data_bolus) - res1_iv_c <- advan("1cmt_iv_bolus", cpp=TRUE)(data_bolus) - res1_iv_ode <- sim(ode = mod_1cmt, regimen = reg_bolus, parameters = parameters, t_obs = t_obs) - - # AUC R - expect_equal(round(res1_iv_r[res1_iv_r$TIME %in% t_obs,]$AUC, 5), round(res1_iv_ode[res1_iv_ode$comp == 2,]$y, 5)) - - #AUC-C - expect_equal(round(res1_iv_c[res1_iv_c$TIME %in% t_obs,]$AUC, 5), round(res1_iv_ode[res1_iv_ode$comp == 2,]$y, 5)) -}) - -test_that("Two compartment bolus ADVAN runs", { - skip_on_cran() - res2_iv_r <- advan("2cmt_iv_bolus", cpp=FALSE)(data_bolus) - res2_iv_c <- advan("2cmt_iv_bolus", cpp=TRUE)(data_bolus) - res2_iv_ode <- sim(ode = mod_2cmt, regimen = reg_bolus, parameters = parameters, t_obs = t_obs) - expect_equal( - round(res2_iv_r[res2_iv_r$TIME %in% t_obs,]$AUC, 5), - round(res2_iv_c[res2_iv_c$TIME %in% t_obs,]$AUC, 5) - ) - # AUC R - expect_equal( - round(res2_iv_r[res2_iv_r$TIME %in% t_obs,]$AUC, 5), - round(res2_iv_ode[res2_iv_ode$comp == 3,]$y, 5) - ) - - #AUC-C - expect_equal( - round(res2_iv_c[res2_iv_c$TIME %in% t_obs,]$AUC, 5), - round(res2_iv_ode[res2_iv_ode$comp == 3,]$y, 5) - ) -}) - -test_that("Two compartment infusion ADVAN runs", { - skip_on_cran() - res2_inf_r <- advan("2cmt_iv_infusion", cpp=FALSE)(data_infusion) - res2_inf_c <- advan("2cmt_iv_infusion", cpp=TRUE)(data_infusion) - res2_inf_ode <- sim(ode = mod_2cmt, regimen = reg_infusion, parameters = parameters, t_obs = t_obs) - - expect_equal( - round(res2_inf_r[res2_inf_r$TIME %in% t_obs,]$AUC, 5), - round(res2_inf_c[res2_inf_c$TIME %in% t_obs,]$AUC, 5) - ) - - # AUC R - expect_equal( - round(res2_inf_r[res2_inf_r$TIME %in% t_obs,]$AUC, 5), - round(res2_inf_ode[res2_inf_ode$comp == 3,]$y, 5) - ) - - #AUC-C - expect_equal( - round(res2_inf_c[res2_inf_c$TIME %in% t_obs,]$AUC, 5), - round(res2_inf_ode[res2_inf_ode$comp == 3,]$y, 5) - ) - -}) - -test_that("Three compartment bolus ADVAN runs", { - skip_on_cran() - res3_iv_r <- advan("3cmt_iv_bolus", cpp=FALSE)(data_bolus) - res3_iv_c <- advan("3cmt_iv_bolus", cpp=TRUE)(data_bolus) - res3_iv_ode <- sim(ode = mod_3cmt, regimen = reg_bolus, parameters = parameters, t_obs = t_obs) - expect_equal( - round(res3_iv_r[res3_iv_r$TIME %in% t_obs,]$AUC, 5), - round(res3_iv_c[res3_iv_c$TIME %in% t_obs,]$AUC, 5) - ) - # AUC R - expect_equal( - round(res3_iv_r[res3_iv_r$TIME %in% t_obs,]$AUC, 5), - round(res3_iv_ode[res3_iv_ode$comp == 4,]$y, 5) - ) - - #AUC-C - expect_equal( - round(res3_iv_c[res3_iv_c$TIME %in% t_obs,]$AUC, 5), - round(res3_iv_ode[res3_iv_ode$comp == 4,]$y, 5) - ) -}) - -test_that("Three compartment iv ADVAN runs", { - skip_on_cran() - res3_iv_r <- advan("3cmt_iv_infusion", cpp=FALSE)(data_infusion) - res3_iv_c <- advan("3cmt_iv_infusion", cpp=TRUE)(data_infusion) - res3_iv_ode <- sim(ode = mod_3cmt, regimen = reg_infusion, parameters = parameters, t_obs = t_obs) - expect_equal( - round(res3_iv_r[res3_iv_r$TIME %in% t_obs,]$AUC, 5), - round(res3_iv_c[res3_iv_c$TIME %in% t_obs,]$AUC, 5) - ) - # AUC R - expect_equal( - round(res3_iv_r[res3_iv_r$TIME %in% t_obs,]$AUC, 5), - round(res3_iv_ode[res3_iv_ode$comp == 4,]$y, 5) - ) - - #AUC-C - expect_equal( - round(res3_iv_c[res3_iv_c$TIME %in% t_obs,]$AUC, 5), - round(res3_iv_ode[res3_iv_ode$comp == 4,]$y, 5) - ) -}) - diff --git a/tests/testthat/test_advan_with_covariates.R b/tests/testthat/test_advan_with_covariates.R deleted file mode 100644 index 7f838358..00000000 --- a/tests/testthat/test_advan_with_covariates.R +++ /dev/null @@ -1,62 +0,0 @@ -test_that("Analytic and ODE models with covariates are the same", { - skip_on_cran() - - ## Create dataset - dose <- 100 - interval <- 12 - n_days <- 2 - parameters <- list( - CL = 10, - V = 50, - KA = 0.5, - Q = 5, - V2 = 100, - Q2 = 3, - V3 = 150, - F1 = 1 - ) - t_obs <- seq(0, 40, .1) - reg_bolus <- new_regimen( - amt = dose, - times = seq(0, interval * n_days * (24/interval), interval), - type = "bolus" - ) - ## there is slight difference in how bolus doses are handled. - ## Analytical equation is perhaps more consistent, so not testing - ## simulations at dose times. Should look into later. - t_obs <- t_obs[! t_obs %in% reg_bolus$dose_times] - covariates <- list(WT = new_covariate(80), CRCL=new_covariate(4.5)) - - ## Using analytic equations model: - data_ana <- sim( - analytical = "1cmt_iv_bolus", - parameters = parameters, - covariates = covariates, - regimen = reg_bolus, - t_obs = t_obs, - covariate_model = "CL = CL * (CRCL / 3)^0.75; V = V * (WT / 70.0)" - ) - - ## Using ODE model: - mod1 <- new_ode_model( - code = " - dAdt[1] = -( (CL*pow(CRCL/3.0, 0.75)) / (V*WT/70.0) ) * A[1]; - ", - covariates = covariates, - obs = list(cmt = 1, scale = "V*WT/70.0"), dose = list(cmt = 1) - ) - data_ode <- sim( - ode = mod1, - parameters = parameters, - covariates = covariates, - regimen = reg_bolus, - t_obs = t_obs, - duplicate_t_obs = TRUE, - only_obs = TRUE - ) - - expect_equal(nrow(data_ana), nrow(data_ode)) - expect_equal(round(data_ana$y,4), round(data_ode$y, 4)) -}) - - diff --git a/tests/testthat/test_bioav_def.R b/tests/testthat/test_bioav_def.R deleted file mode 100644 index c4d7674a..00000000 --- a/tests/testthat/test_bioav_def.R +++ /dev/null @@ -1,44 +0,0 @@ -parameters <- list(KA = 0.5, CL = 5, V = 50) -covs <- list(WT = new_covariate(50)) -reg <- new_regimen(amt = 100, n = 1, interval = 12, type = "bolus") -mod1 <- oral_1cmt_allometric # defined in setup.R using same covs and pars as above -y1 <- sim_ode(mod1, parameters = parameters, regimen = reg, covariates = covs, only_obs=TRUE)$y - -test_that("bioav numeric option working", { - skip_on_cran() - mod2 <- new_ode_model( - code = " - CLi = CL * pow(WT/70, 0.75) - dAdt[1] = -KA * A[1] - dAdt[2] = KA*A[1] - (CLi/V)*A[2] - ", - dose = list(cmt = 1, bioav = 0.5), - covariates = covs, - declare_variables = "CLi", - parameters = parameters - ) - - y2 <- sim_ode(mod2, parameters = parameters, regimen = reg, covariates = covs, only_obs=TRUE)$y - - expect_equal(round(y1,1), round(y2*2, 1)) -}) - -test_that("bioav string option working", { - skip_on_cran() - mod3 <- new_ode_model( - code = " - CLi = CL * pow(WT/70, 0.75) - dAdt[1] = -KA * A[1] - dAdt[2] = KA*A[1] - (CLi/V)*A[2] - ", - dose = list(cmt = 1, bioav = "WT/70"), - covariates = covs, - declare_variables = "CLi", - parameters = parameters - ) - - y3 <- sim_ode(mod3, parameters = parameters, regimen = reg, covariates = covs, only_obs=TRUE)$y - - expect_equal(round(y1*(50/70),1), round(y3, 1)) -}) - diff --git a/tests/testthat/test_calc_ss_analytic.R b/tests/testthat/test_calc_ss_analytic.R index 43634725..54ca63c5 100644 --- a/tests/testthat/test_calc_ss_analytic.R +++ b/tests/testthat/test_calc_ss_analytic.R @@ -1,3 +1,7 @@ +# Uses models defined in setup.R: +# - mod_1cmt_iv, mod_2cmt_iv, mod_1cmt_oral +# - pk_3cmt_iv, pk_2cmt_oral, pk_3cmt_oral (conditional, NOT_CRAN only) + # shared parameters dose <- 100 interval <- 12 @@ -6,13 +10,6 @@ reg_oral <- new_regimen(amt = dose, interval = interval, n = n_ss, type = "oral" reg_bolus <- new_regimen(amt = dose, interval = interval, n = n_ss, type = "bolus") reg_inf <- new_regimen(amt = dose, interval = interval, n = n_ss, type = "infusion") t_obs <- max(reg_oral$dose_times) + interval -# Uses models defined in setup.R - -if (identical(Sys.getenv("NOT_CRAN"), "true")) { - pk_3cmt_iv <- new_ode_model("pk_3cmt_iv") -} - -#delta <- function(x, ref) { abs(x-ref)/ref } test_that("1 cmt oral", { par <- list(CL = 5, V = 100, KA = 1) @@ -72,7 +69,6 @@ test_that("1-cmt iv infusion", { test_that("2-cmt oral", { skip_on_cran() - pk_2cmt_oral <- new_ode_model("pk_2cmt_oral") par <- list(CL = 5, V = 100, Q = 3, V2 = 150, KA = 1) res_ana <- calc_ss_analytic(f = "2cmt_oral", dose = dose, interval = interval, parameters = par) res_ode <- sim(pk_2cmt_oral, parameters = par, regimen = reg_oral, t_obs = t_obs, only_obs = F)$y @@ -98,7 +94,6 @@ test_that("2-cmt infusion", { test_that("3-cmt oral", { skip_on_cran() - pk_3cmt_oral <- new_ode_model("pk_3cmt_oral") par <- list(CL = 5, V = 100, Q = 3, V2 = 150, Q2 = 6, V3 = 250, KA = 1) res_ana <- calc_ss_analytic(f = "3cmt_oral", dose = dose, interval = interval, parameters = par) res_ode <- sim(pk_3cmt_oral, parameters = par, regimen = reg_oral, t_obs = t_obs, only_obs = F)$y diff --git a/tests/testthat/test_calc_parameters.R b/tests/testthat/test_calculate_parameters.R similarity index 100% rename from tests/testthat/test_calc_parameters.R rename to tests/testthat/test_calculate_parameters.R diff --git a/tests/testthat/test_check_mixture_model.R b/tests/testthat/test_check_mixture_model.R new file mode 100644 index 00000000..8e97be4f --- /dev/null +++ b/tests/testthat/test_check_mixture_model.R @@ -0,0 +1,17 @@ +test_that("expected mixture model errors are caught", { + mixture1 <- list(CL = list(values = c(5, 10), probability = 0.3)) + mixture2 <- list( + CL = list(values = c(5, 10), probability = 0.3), + V = list(values = c(10, 20), probability = 0.5) + ) + mixture_multi <- list( + CL = list(values = c(5, 6, 7), probability =c(0.3, 0.3, 0.4)) + ) + mixture_oops1 <- list(CL = list(values = c(5, 10), probability = 30)) + mixture_oops2 <- list(CL = list(values = c(5, 10), probability = -0.3)) + expect_error(check_mixture_model(mixture2, c("CL", "V"))) + expect_error(check_mixture_model(mixture1, c("Ke", "V"))) + expect_error(check_mixture_model(mixture_multi, c("CL", "V"))) + expect_error(check_mixture_model(mixture_oops1, c("CL", "V"))) + expect_error(check_mixture_model(mixture_oops2, c("CL", "V"))) +}) diff --git a/tests/testthat/test_cmt_mapping.R b/tests/testthat/test_cmt_mapping.R deleted file mode 100644 index 48af64e6..00000000 --- a/tests/testthat/test_cmt_mapping.R +++ /dev/null @@ -1,56 +0,0 @@ -pk1cmt_oral_code <- new_ode_model( - code = "dAdt[1] = -KA*A[1]; dAdt[2] = KA*A[1] - (CL/V)*A[2];", - obs = list(cmt = 2, scale="V"), - cmt_mapping = list(oral = 1, infusion = 2, bolus = 2) -) - -test_that("Compartment mapping is added to attributes", { - expect_equal(attr(pk1cmt_oral_code, "cmt_mapping")[["oral"]], 1) - expect_equal(attr(pk1cmt_oral_code, "cmt_mapping")[["infusion"]], 2) -}) - -test_that("Admin route is interpreted and simulated correctly", { - regimen <- new_regimen( - amt = c(100, 100, 100, 100), - times = c(0, 12, 24, 36), - type = c("oral", "oral", "infusion", "infusion"), - t_inf = c(0, 0, 1, 1) - ) - p <- list(KA = 1, CL = 5, V = 50) - res <- sim_ode( - ode = pk1cmt_oral_code, - parameters = p, - regimen = regimen, - only_obs = FALSE - ) - expect_equal(round(res$y[res$comp == 1 & res$t == 25], 4), 2e-04) - expect_true(res$y[res$comp == 2 & res$t == 25] >= 100) -}) - -test_that("multiple scaling types on one compartment works", { - skip_on_cran() - mod <- new_ode_model( - code = " - dAdt[1] = -KA * A[1]; - dAdt[2] = -(CL/V) * A[2] + KA*A[1]; - ", - obs = list( - cmt = c(2, 2), - scale = c(1, "V"), - label = c("abs", "conc") - ), - cpp_show_code = FALSE - ) - par <- list(CL = 5, V = 50, KA = .5) - reg <- new_regimen(amt = 100, n = 5, interval = 12) - res <- sim_ode( - ode = mod, - parameters = par, - regimen = reg, - only_obs = TRUE - ) - dat <- cbind(res[res$comp == "abs",]$y, res[res$comp == "conc",]$y) - expect_true("PKPDsim" %in% class(mod)) - expect_equal(length(unique(res$comp)), 2) - expect_equal(round(dat[,1],1), round(dat[,2]*par$V,1)) -}) diff --git a/tests/testthat/test_compare_results.R b/tests/testthat/test_compare_results.R deleted file mode 100644 index 73652e9f..00000000 --- a/tests/testthat/test_compare_results.R +++ /dev/null @@ -1,266 +0,0 @@ -## models: shared between tests and take a while to compile -# - oral models -## Uses model defined in setup.R -pk1cmt_oral_anal = function(t, dose, KA, V, CL) { - dose*KA/(V*(KA-CL/V))*(exp(-(CL/V) * t)-exp(-KA * t)) -} -pk1cmt_oral_code <- new_ode_model( - code = "dAdt[1] = -KA*A[1]; dAdt[2] = KA*A[1] - (CL/V)*A[2];", - obs=list(cmt = 2, scale="V") -) - -# - iv models -## Uses model defined in setup.R - -# - model with dose cmt specified -dose_in_cmt_2 <- new_ode_model( - code = " - dAdt[1] = -KA * A[1]; - dAdt[2] = KA*A[1] -(CL/V) * A[2] - dAdt[3] = S2*(A[2]-A[3]) - ", - obs = list(cmt=2, scale="V"), - dose = list(cmt = 2), - cpp_show_code = FALSE -) - - -test_that("Library and custom C++ and code matches analytic soln", { - p <- list(KA = 1, CL = 5, V = 50) - t_obs <- c(0:72) - t_obs2 <- t_obs + 0.1234 # also needs to be producing results with non-integer times - dose <- 100 - t_dose <- c(0) - regimen <- new_regimen(amt=dose, times = t_dose, type = "oral") - - pk1cmt_oral_lib <- sim_ode( - ode = mod_1cmt_oral, - parameters = p, - regimen = regimen, - t_obs = t_obs, - int_step_size = 0.1, - duplicate_t_obs = TRUE, - only_obs=TRUE - ) - - pk1cmt_oral_code <- sim_ode( - ode = pk1cmt_oral_code, - parameters = p, - duplicate_t_obs = TRUE, - regimen=regimen, - t_obs=t_obs, - int_step_size = 0.1, - only_obs=TRUE - ) - - pk1cmt_oral_anal_res <- pk1cmt_oral_anal(t_obs, dose, p$KA, p$V, p$CL) - expect_equal(round(pk1cmt_oral_lib$y, 3), round(pk1cmt_oral_anal_res, 3)) - expect_equal(round(pk1cmt_oral_code$y, 3), round(pk1cmt_oral_anal_res, 3)) -}) - - -test_that("precision in time does not impact # obs returned", { - regimen_mult <- new_regimen( - amt = rep(12.8, 3), - times = c(0, 6, 12), - type="infusion", - t_inf = 2 - ) - t_obs <- c(11.916, 14.000, 16.000, 17.000, 30) - tmp <- sim_ode( - ode = mod_1cmt_iv, - parameters = list(CL = 5, V = 50), - regimen = regimen_mult, - t_obs = t_obs, - only_obs = TRUE - ) - expect_equal(tmp$t, t_obs) -}) - -test_that("test bug EmCo 20150925", { - xtim <- c(0, 2, 4, 8, 12, 24) - sujdos <- 320 - param <- list(KA = 1.8, V = 30, CL = 1.7) - regim <- new_regimen(amt = sujdos, times = c(0, 12), type= "bolus") - out <- sim_ode(ode = mod_1cmt_oral, parameters=param, regimen=regim, t_obs = xtim, only_obs = TRUE) - expect_equal(out$t, xtim) -}) - -test_that("model size is appropriate (bug: JeHi 20151204)", { - pk3cmt <- new_ode_model( - code = " - dAdt[1] = -KA*A[1]; - dAdt[2] = KA*A[1] -(Q/V)*A[2] + (Q/V2)*A[3] -(CL/V)*A[2]; - dAdt[3] = -(Q/V2)*A[3] + (Q/V)*A[2]; - ", - obs = list(cmt = 2, scale = "V") - ) - expect_equal( attr(pk3cmt, "size"), 3) -}) - -test_that("Dose is added to correct compartment: specified by model", { - set.seed(90) - p <- list(CL = 1, V = 10, KA = 0.5, S2 = .1) - r <- new_regimen(amt = 100, times = c(0), type = "infusion") - dat <- sim_ode( - ode = dose_in_cmt_2, - n_ind = 1, - omega = cv_to_omega(par_cv = list("CL"=0.1, "V"=0.1, "KA" = .1), p), - parameters = p, - regimen = r, - verbose = FALSE, - t_max = 48 - ) - # Dose should be in cmt 2 - expect_equal(dat$y[dat$comp == 1], rep(0, 50)) - expect_true(all(dat$y[dat$comp == 2][-1] > 0)) -}) - -test_that("Dose is added to correct compartment: override model by regimen", { - set.seed(60) - p <- list(CL = 1, V = 10, KA = 0.5, S2 = .1) - r <- new_regimen( - amt = c(100, 100, 100), - times = c(0, 6, 12), - cmt = c(1,2,3), - type = "bolus" - ) - dat <- sim_ode( - ode = dose_in_cmt_2, - n_ind = 1, - omega = cv_to_omega(par_cv = list("CL"=0.1, "V"=0.1, "KA" = .1), p), - parameters = p, - regimen = r, - verbose = FALSE, - t_max = 48 - ) - # Dose should be in cmt 1, 2 and 3 - expect_true(all(dat$y[dat$comp == 1 & dat$t > 0] > 0)) - expect_true(max(diff(dat$y[dat$comp == 2])) > 95) - expect_true(max(diff(dat$y[dat$comp == 3])) > 95) -}) - -test_that("Infusion works for all compartments", { - set.seed(44) - # Part 1: Specify cmt only with model - p <- list(CL = 1, V = 10, KA = 0.5, S2 = .1) - r <- new_regimen( - amt = c(100, 100, 100), - times = c(0, 6, 12), - cmt = c(1,2,3), - t_inf = 3, - type = "infusion" - ) - dat <- sim_ode( - ode = dose_in_cmt_2, - n_ind = 1, - omega = cv_to_omega(par_cv = list("CL"=0.1, "V"=0.1, "KA" = .1), p), - parameters = p, - regimen = r, - verbose = FALSE, - t_max = 48 - ) - expect_true(all(dat$y[dat$comp == 1 & dat$t > 0 ] > 0)) - expect_true(max(diff(dat$y[dat$comp == 2])) > 25) - expect_true(max(diff(dat$y[dat$comp == 3])) > 25) - expect_equal(round(max(dat$y[dat$comp == 2]), 1), 131.2) - expect_equal(round(max(dat$y[dat$comp == 3]), 1), 148.4) -}) - -test_that("Duplicate obs returned when specified in arg", { - # for first 2 doses, infusion time will just be ignored, but a value has to be specified in the vector - p <- list(CL = 1, V = 10, KA = 0.5, S2=.1) - r <- new_regimen( - amt = c(100, 100, 100, 100), - times = c(0, 6, 12, 18), - cmt = c(2, 2, 1, 1), - t_inf = c(1, 1, 1, 1), - type = c("bolus", "bolus", "infusion", "infusion") - ) - dat <- sim_ode( - ode = mod_1cmt_oral, - n_ind = 1, - omega = cv_to_omega(par_cv = list("CL"=0.1, "V"=0.1, "KA" = .1), p), - parameters = p, - regimen = r, - t_obs = c(1, 2, 3, 4, 4, 4, 6), ## see duplicate obs here - duplicate_t_obs = T, - only_obs = FALSE - ) - expect_equal(length(dat[dat$t == 4,]$y), 9) - expect_equal(length(dat$y), 21) - expect_equal(sum(is.na(dat$y)), 0) -}) - -test_that("Custom t_obs is returned", { - t_obs <- seq(from = 0, to = 24, by = .1) - p <- list(CL = 1, V = 10, KA = 0.5, S2=.1) - r <- new_regimen( - amt = c(100, 100, 100, 100), - times = c(0, 6, 12, 18), - cmt = c(2, 2, 1, 1), - t_inf = c(1, 1, 1, 1), - type = c("bolus", "bolus", "infusion", "infusion") - ) - dat <- sim_ode( - ode = mod_1cmt_oral, - n_ind = 1, - omega = cv_to_omega(par_cv = list("CL"=0.1, "V"=0.1, "KA" = .1), p), - parameters = p, - regimen = r, - t_obs = t_obs, - only_obs = T - ) - expect_equal(mean(diff(t_obs)), mean(diff(dat$t))) -}) - -test_that("if covariate time is at end of infusion, end of infusion is still recorded", { - # Bug reported by JF - pop_est <- list(CL = 1.08, V = 0.98) - regimen <- new_regimen( - amt = c(1500, 1000, 1500, 1500, 1500, 1500, 1500), - type = "infusion", - t_inf = c(2, 1, 2, 2, 1, 1, 1), - times = c(0, 10.8666666666667, 20.4333333333333, 32.0666666666667, 46.9, 54.9, 62.9 ) - ) - covs <- list( - WT = new_covariate(value = c(60, 65), times = c(0, 47.9)), - CRCL = new_covariate(8), CVVH = new_covariate(0) - ) - pksim <- sim( - ode = mod_1cmt_iv, - parameters = pop_est, - covariates = covs, - regimen = regimen, - checks = TRUE, - only_obs = TRUE - ) - expect_true(all(pksim$y < 1000)) -}) - -test_that("Covariate table simulation runs", { - # this test used to be in the covariate_table_to_list file but - # makes more sense here. - p <- list(CL = 5, V = 50) - reg <- new_regimen (amt = 100, n = 4, interval = 12, type = "bolus", cmt=1) - om <- c(0.01, 1, 0.01) - cov_table <- data.frame( - id=c(1, 1, 2, 3), - WT = c(40, 45, 50, 60), - SCR = c(50, 150, 90,110), - t = c(0, 5, 0, 0) - ) - - dat <- sim( - mod_1cmt_iv, - parameters = p, - regimen = reg, - covariates_table = cov_table, - covariates_implementation = list(SCR = "interpolate"), - omega = NULL, - n_ind = 3, - only_obs = T, - output_include = list(parameters = TRUE, covariates=TRUE) - ) - expect_equal(length(unique(dat$id)), 3) -}) diff --git a/tests/testthat/test_get_var_y.R b/tests/testthat/test_get_var_y.R index 264843f6..750a5592 100644 --- a/tests/testthat/test_get_var_y.R +++ b/tests/testthat/test_get_var_y.R @@ -1,6 +1,9 @@ +# Uses models defined in setup.R: +# - mod_1cmt_iv +# - mod_2cmt_iv +# - mod_1cmt_iv_mm (conditional, NOT_CRAN only) ## Set up simulations to test variance: -## Uses model defined in setup.R reg <- new_regimen( amt = 100, n = 3, @@ -12,7 +15,7 @@ reg <- new_regimen( par <- list(CL = 5, V = 50) omega <- c(0.1, 0.0, 0.1) t_obs <- c(2, 48) -res <- sim_ode( +res <- sim( mod_1cmt_iv, parameters = par, regimen = reg, @@ -29,7 +32,8 @@ test_that("delta approximation and full simulation match", { t_obs = t_obs, regimen = reg, omega = omega, - method = "delta" + method = "delta", + seed = 12345 ) v_sim <- get_var_y( model = mod_1cmt_iv, @@ -39,7 +43,8 @@ test_that("delta approximation and full simulation match", { regimen = reg, omega = omega, method="sim", - n_ind = 2500 + n_ind = 2500, + seed = 12345 ) expect_true(max(abs(v_delta$regular - v_sim$regular)/v_sim$regular) <0.1) expect_true(max(abs(v_delta$log - v_sim$log)/v_sim$log) < 0.15) @@ -63,7 +68,8 @@ test_that("Confidence interval instead of SD", { regimen = reg2, omega = omega, method = "delta", - q = CI_range + q = CI_range, + seed = 12345 ) v2_sim <- get_var_y( model = mod_1cmt_iv, @@ -73,7 +79,8 @@ test_that("Confidence interval instead of SD", { omega = omega, method ="sim", n_ind = 2500, - q = CI_range + q = CI_range, + seed = 12345 ) expect_equal(dim(v1_delta$regular), c(2, 2)) expect_equal(dim(v2_sim$regular), c(2, 2)) @@ -88,7 +95,6 @@ test_that("Confidence interval instead of SD", { expect_true(max(abs(v1_delta_v - v2_sim_v )/v2_sim_v) < 0.1) }) - test_that("Two compartment model", { skip_on_cran() set.seed(80) @@ -111,7 +117,8 @@ test_that("Two compartment model", { parameters = par2, t_obs = t_obs, regimen = reg, - omega = omega2 + omega = omega2, + seed = 12345 ) v2 <- get_var_y( model = mod_2cmt_iv, @@ -120,7 +127,8 @@ test_that("Two compartment model", { regimen = reg, omega = omega2, method = "sim", - n_ind = 2000 + n_ind = 2000, + seed = 12345 ) expect_true(all(abs(((v1$regular - v2$regular)/v2$regular)) < 0.05)) @@ -130,13 +138,12 @@ test_that("Two compartment model", { test_that("One compartment with MM kinetics", { skip_on_cran() - mod3 <- new_ode_model("pk_1cmt_iv_mm") par3 <- list(VMAX = 5, KM = 5, V = 10) omega3 <- c(0.1, 0.05, 0.1, 0.01, 0.01, 0.1) res <- sim_ode( - mod3, + mod_1cmt_iv_mm, parameters = par3, t_obs = t_obs, regimen = reg, @@ -144,20 +151,22 @@ test_that("One compartment with MM kinetics", { ) v1 <- get_var_y( - model = mod3, + model = mod_1cmt_iv_mm, parameters = par3, t_obs = t_obs, regimen = reg, - omega = omega3 + omega = omega3, + seed = 12345 ) v2 <- get_var_y( - model = mod3, + model = mod_1cmt_iv_mm, parameters = par3, t_obs = t_obs, regimen = reg, omega = omega3, method = "sim", - n_ind = 2000 + n_ind = 2000, + seed = 12345 ) expect_true(all(abs(v1$regular - v2$regular)/res$y < 0.05)) diff --git a/tests/testthat/test_iov.R b/tests/testthat/test_iov.R deleted file mode 100644 index 650d1388..00000000 --- a/tests/testthat/test_iov.R +++ /dev/null @@ -1,259 +0,0 @@ -pars <- list( - "kappa_CL_1" = 0, - "kappa_CL_2" = 0, - "kappa_CL_3" = 0, - "eta_CL" = 0, - "CL" = 5, - "V" = 50, - "KA" = 1 -) -pars0 <- list( - "CL" = 5, - "V" = 50, - "KA" = 1 -) -pk0 <- new_ode_model( # no IOV - code = " - dAdt[1] = -KA * A[1] - dAdt[2] = +KA * A[1] -(CL/V) * A[2] - ", - obs = list(cmt = 2, scale = "V"), - dose = list(cmt = 1, bioav = 1), - parameters = names(pars0), - cpp_show_code = F -) -pk1 <- new_ode_model( - code = " - CL_iov = CL * exp(kappa_CL + eta_CL); - dAdt[1] = -KA * A[1] - dAdt[2] = +KA * A[1] -(CL_iov/V) * A[2] - ", - iov = list( - cv = list(CL = 0.2), - n_bins = 3 - ), - obs = list(cmt = 2, scale = "V"), - dose = list(cmt = 1, bioav = 1), - declare_variables = c("kappa_CL", "CL_iov"), - parameters = names(pars), - cpp_show_code = F -) -reg1 <- new_regimen( - amt = 100, - interval = 24, - n = 5, - type = "infusion" -) -iov_var <- 0.3 ^ 2 # 30% IOV - -test_that("Throws error when `iov_bins` supplied but not present in model", { - expect_error({ - sim( - ode = pk0, - parameters = pars0, - regimen = reg1, - omega = c( - 0.3 # IIV in CL - ), - n = 1, - iov_bins = c(0, 24, 48, 72, 999), # !! - omega_type = "normal", - only_obs = TRUE, - output_include = list(parameters = TRUE, variables = TRUE) - ) - }, "No IOV implemented for this model") -}) - -test_that("Throws error when number of `iov_bins` is higher than allowed for model", { - expect_error({ - sim( - ode = pk1, - parameters = pars, - regimen = reg1, - omega = c( - iov_var, # IOV in CL - 0, iov_var, - 0, 0, iov_var, - 0, 0, 0, iov_var, - 0, 0, 0, 0, 0.3 # IIV in CL - ), - n = 1, - iov_bins = c(0, 24, 48, 72, 999), # one bin too many - omega_type = "normal", - only_obs = TRUE, - output_include = list(parameters = TRUE, variables = TRUE) - ) - }, "Number of allowed IOV bins for model is lower") -}) - -test_that("Throws warning when number of `iov_bins` is lower than allowed for model", { - expect_warning({ - sim( - ode = pk1, - parameters = pars, - regimen = reg1, - omega = c( - iov_var, # IOV in CL - 0, iov_var, - 0, 0, iov_var, - 0, 0, 0, iov_var, - 0, 0, 0, 0, 0.3 # IIV in CL - ), - n = 1, - iov_bins = c(0, 24, 999), # one bin too few - omega_type = "normal", - only_obs = TRUE, - output_include = list(parameters = TRUE, variables = TRUE) - )}, "Number of allowed IOV bins for model is higher" - ) -}) - -test_that("IOV is added to parameters", { - skip_on_cran() - set.seed(32) - - dat <- sim( - ode = pk1, - parameters = pars, - regimen = reg1, - omega = c( - iov_var, # IOV in CL - 0, iov_var, - 0, 0, iov_var, - 0, 0, 0, iov_var, - 0, 0, 0, 0, 0.3 # IIV in CL - ), - n = 5, - iov_bins = c(0, 24, 48, 72), - omega_type = "normal", - only_obs = TRUE, - output_include = list(parameters = TRUE, variables = TRUE) - ) - - expect_equal( - signif(dat$kappa_CL[dat$t == 12], 4), - signif(dat$kappa_CL_1[dat$t == 1], 4) - ) - expect_equal( - signif(dat$kappa_CL[dat$t == 36], 4), - signif(dat$kappa_CL_2[dat$t == 1], 4) - ) - expect_equal( - signif(dat$kappa_CL[dat$t == 48], 4), - signif(dat$kappa_CL_3[dat$t == 1], 4) - ) - expect_equal( - signif(exp(dat$kappa_CL + dat$eta_CL) * dat$CL, 4), - signif(dat$CL_iov, 4) - ) -}) - -test_that("Change in F in 2nd bin is applied in 2nd bin and not later.", { - # Previously this was an issue because F, when defined in pk_code(), was not updated before the - # dose was applied to the state vector, so the bioavailability was not applied at the right time. - # This was fixed by rearranging the order of execution in sim.cpp in the main loop. - - pars <- list( - "CL" = 5, - "V" = 50, - "KA" = 1, - "F" = 1 - ) - pk1 <- new_ode_model( - code = " - dAdt[1] = -KA * A[1] - dAdt[2] = +KA * A[1] -(CLi/V) * A[2] - ", - pk_code = " - Fi = F * exp(kappa_F); - CLi = CL; - ", - iov = list( - cv = list(F = 0.2), - n_bins = 3 - ), - obs = list(cmt = 2, scale = "V"), - dose = list(cmt = 1, bioav = "Fi"), - declare_variables = c("kappa_F", "Fi", "CLi"), - parameters = names(pars), - cpp_show_code = F - ) - reg <- new_regimen(amt = 800, interval = 24, n = 10, type = "oral") - - # For a first simulation, we're simulating with no variability across the IOV bins: - pars$kappa_F_1 <- 0 - pars$kappa_F_2 <- 0 - pars$kappa_F_3 <- 0 - args_sim1 <- args <- list( - ode = pk1, - parameters = pars, - regimen = reg, - only_obs = TRUE, - t_obs = seq(0, 50, .25), - iov_bins = c(0L, 24L, 48L, 9999L) - ) - # For a second simulation, we're applying a change in parameter for the 2nd bin (24-48 hrs). - # This should affect predictions from 24 onward. - pars$kappa_F_2 <- 1 # 2nd bin - args_sim2 <- args <- list( - ode = pk1, - parameters = pars, - regimen = reg, - only_obs = TRUE, - t_obs = seq(0, 50, .25), - iov_bins = c(0L, 24L, 48L, 9999L) - ) - res1 <- do.call("sim_ode", args = args_sim1) - res2 <- do.call("sim_ode", args = args_sim2) - expect_true(min(res1[res1$y != res2$y,]$t) <= 25) -}) - -test_that("error is not invoked when using parameters_table", { - parameters_table <- data.frame( - CL = runif(10), V = runif(10), KA = runif(10), eta_CL = runif(10), - kappa_CL_1 = 0, kappa_CL_2 = 0, kappa_CL_3 = 0, kappa_CL_4 = 0 - ) - - # specifying both parameters_table but for a model with IOV should not fail! - expect_silent( - dat <- sim( - ode = pk1, - parameters_table = parameters_table, - regimen = reg1, - omega = c( - iov_var, # IOV in CL - 0, iov_var, - 0, 0, iov_var, - 0, 0, 0, iov_var, - 0, 0, 0, 0, 0.3 # IIV in CL - ), - n = 5, - iov_bins = c(0, 24, 48, 72), - omega_type = "normal", - only_obs = TRUE, - output_include = list(parameters = TRUE, variables = TRUE) - ) - ) - - # specifying both parameters and parameters_table should fail - expect_error( - dat <- sim( - ode = pk1, - parameters = pars, - parameters_table = parameters_table, - regimen = reg1, - omega = c( - iov_var, # IOV in CL - 0, iov_var, - 0, 0, iov_var, - 0, 0, 0, iov_var, - 0, 0, 0, 0, 0.3 # IIV in CL - ), - n = 5, - iov_bins = c(0, 24, 48, 72), - omega_type = "normal", - only_obs = TRUE, - output_include = list(parameters = TRUE, variables = TRUE) - ) - ) -}) diff --git a/tests/testthat/test_mixture_model.R b/tests/testthat/test_mixture_model.R deleted file mode 100644 index 4289592a..00000000 --- a/tests/testthat/test_mixture_model.R +++ /dev/null @@ -1,64 +0,0 @@ -## Setup model + params -covs <- list(WT = PKPDsim::new_covariate(70)) -mod <- new_ode_model( - code = " - dAdt[0] = -(CL*(WT/70.0)/V)*A[0]; - ", - pk_code = " ", - obs = list(cmt = 1, scale = "V"), - mixture = list(CL = list(values = c(5, 15), probability = 0.3)), - covariates = covs, -) -par <- list(CL = 3, V = 50) -reg <- new_regimen(amt = 250, n = 5, interval = 6, type = 'infusion', t_inf = 1) -t_obs <- seq(0, 36, 4) - -test_that("mixture model works properly for single patient", { - res0 <- sim_ode(mod, parameters = par, regimen = reg, covariates = covs, t_obs = t_obs, only_obs=T) # mixture_group not supplied - res1 <- sim(mod, parameters = par, regimen = reg, t_obs = t_obs, covariates = covs, mixture_group = 1, only_obs=T) - res2 <- sim(mod, parameters = par, regimen = reg, t_obs = t_obs, covariates = covs, mixture_group = 2, only_obs=T) - expect_equal(round(res0[res0$t == 24,]$y, 2), 9.07) # should use whatever is in `parameters` - expect_equal(round(res1[res1$t == 24,]$y, 2), 5.82) - expect_equal(round(res2[res2$t == 24,]$y, 2), 1.15) -}) - -test_that("mixture model works properly when vectorized (using parameters_table)", { - partab <- data.frame(CL = rep(0, 6), V = rep(50, 6)) - suppressMessages({ - expect_error(sim_ode(mod, parameters_table = partab, regimen = reg, t_obs = t_obs, covariates = covs, mixture_group = 1, only_obs=T)) - res1 <- sim(mod, parameters_table = partab, regimen = reg, t_obs = t_obs, covariates = covs, mixture_group = rep(1, 6), only_obs=T) - res2 <- sim(mod, parameters_table = partab, regimen = reg, t_obs = t_obs, covariates = covs, mixture_group = rep(c(1,2), 3), only_obs=T, output_include = list(parameters = TRUE)) - }) - expect_equal(round(res1[res1$t == 24,]$y, 2), rep(5.82, 6)) - expect_equal(round(res2[res2$t == 24,]$y, 2), rep(c(5.82, 1.15), 3)) - expect_equal(res2[res2$id == 1,]$CL[1], 5) - expect_equal(res2[res2$id == 2,]$CL[1], 15) - expect_equal(res2[res2$id == 3,]$CL[1], 5) -}) - -test_that("mixture model works properly when vectorized (using covariates_table)", { - covtab <- data.frame(ID = 1:8, WT = rep(seq(40, 130, 30), 2)) - suppressMessages({ - expect_error(sim(mod, parameters = par, covariates_table = covtab, regimen = reg, t_obs = t_obs, mixture_group = 1, only_obs=T)) - res <- sim(mod, parameters = par, covariates_table = covtab, regimen = reg, t_obs = t_obs, mixture_group = rep(c(1, 2), each=4), only_obs=T) - }) - expect_equal(round(res[res$t == 24,]$y, 2), c(9.39, 5.82, 3.83, 2.65, 2.99, 1.15, 0.52, 0.25)) -}) - -test_that("expected mixture model errors are caught", { - mixture1 <- list(CL = list(values = c(5, 10), probability = 0.3)) - mixture2 <- list( - CL = list(values = c(5, 10), probability = 0.3), - V = list(values = c(10, 20), probability = 0.5) - ) - mixture_multi <- list( - CL = list(values = c(5, 6, 7), probability =c(0.3, 0.3, 0.4)) - ) - mixture_oops1 <- list(CL = list(values = c(5, 10), probability = 30)) - mixture_oops2 <- list(CL = list(values = c(5, 10), probability = -0.3)) - expect_error(check_mixture_model(mixture2, c("CL", "V"))) - expect_error(check_mixture_model(mixture1, c("Ke", "V"))) - expect_error(check_mixture_model(mixture_multi, c("CL", "V"))) - expect_error(check_mixture_model(mixture_oops1, c("CL", "V"))) - expect_error(check_mixture_model(mixture_oops2, c("CL", "V"))) -}) diff --git a/tests/testthat/test_multi_obs.R b/tests/testthat/test_multi_obs.R deleted file mode 100644 index f64511f7..00000000 --- a/tests/testthat/test_multi_obs.R +++ /dev/null @@ -1,236 +0,0 @@ -# Example multiple observation types (e.g. different compartments! not just different residual errors) - -## define parameters -vars <- c("CONC", "METAB", "METAB2", "ACT") -pk1 <- new_ode_model( - code = "dAdt[1] = -(CL/V)*A[1]; CONC = 1000*A[1]/V; METAB = CONC/2; METAB2 = CONC * t; ACT = 15", - obs = list(variable = vars), - declare_variables = vars, - cpp_show_code = F -) - -regimen <- new_regimen(amt = 100, interval = 12, n = 5, type="infusion", t_inf = 1) -parameters <- list("CL" = 15, "V" = 150) -omega <- PKPDsim::cv_to_omega(list("CL" = 0.2, "V" = 0.2), parameters[1:2]) - -test_that("obs types are output by `sim`", { - obs_type <- c(1,2,1,3,1) - data <- sim( - ode = pk1, - parameters = list(CL = 20, V = 200), - regimen = regimen, - int_step_size = 0.1, - only_obs = TRUE, - obs_type = obs_type, - t_obs = c(2, 4, 6, 8, 12), - output_include = list("variables" = TRUE) - ) - expect_equal(data$obs_type, obs_type) - expect_equal(data$y, diag(as.matrix(data[1:5,5+obs_type]))) -}) - - -test_that("check obs at same timepoint but with different obs_type", { - obs_type <- c(1,2,1,3,1) - t_same <- sim( - ode = pk1, - parameters = list(CL = 20, V = 200), - regimen = regimen, - int_step_size = 0.1, - only_obs = TRUE, - obs_type = obs_type, - t_obs = c(2, 4, 4, 8, 8), - output_include = list("variables" = TRUE) - ) - expect_equal(t_same$t[2], t_same$t[3]) - expect_equal(t_same$t[4], t_same$t[5]) - expect_equal(t_same$y[3], t_same$METAB[3]) - expect_equal(t_same$y[2], t_same$CONC[2]) - expect_equal(t_same$y[5], t_same$METAB2[5]) - expect_equal(t_same$y[4], t_same$CONC[4]) -}) - - -test_that("check that residual error correctly applied to right var", { - obs_type <- c(1,2,1,3,1) - set.seed(12345) - ruv_term3 <- list(prop = c(0, 0, 0.1), add = c(0, 0, 0.1)) - ruv_term1 <- list(prop = c(.1, 0, 0), add = c(1, 0, 0)) - - data2 <- sim( - ode = pk1, - parameters = list(CL = 20, V = 200), - regimen = regimen, - int_step_size = 0.1, - only_obs = TRUE, - obs_type = obs_type, - t_obs = c(2, 4, 6, 8, 12), - output_include = list("variables" = TRUE), - res_var = ruv_term3 - ) - y <- diag(as.matrix(data2[1:5, 5 + obs_type])) - - # no residual error for obs_type 1 and 2, only 3 - expect_equal(data2$y[-4], y[-4]) - expect_true(data2$y[-4][4] != y[4]) - - data3 <- sim( - ode = pk1, - parameters = list(CL = 20, V = 200), - regimen = regimen, - int_step_size = 0.1, - only_obs = TRUE, - obs_type = obs_type, - t_obs = c(2, 4, 6, 8, 12), - output_include = list("variables" = TRUE), - res_var = ruv_term1 - ) - y <- diag(as.matrix(data2[1:5, 5 + obs_type])) - - # no residual error for obs_type 2/3, only and 1 - expect_equal(data3$y[-c(1,3,5)], y[-c(1,3,5)]) - expect_true(all(data3$y[c(1,3,5)] != y[c(1,3,5)])) -}) - - -test_that("specifying ruv as multi-type when only 1 obs_type", { - obs_type <- c(1,2,1,3,1) - ruv_single <- list(prop = 0.1, add = 1) - ruv_multi <- list(prop = c(0.1, 1), add = c(1, 20)) - - s_single <- sim( - ode = pk1, - parameters = parameters, - n_ind = 1, - regimen = regimen, - only_obs = TRUE, - t_obs = c(2,4,6,8), - seed = 123456, - res_var = ruv_single - ) - - ## specified as multi, but obs_type is all 1, so should - ## give same results as s_single - s_multi1 <- sim( - ode = pk1, - parameters = parameters, - n_ind = 1, - regimen = regimen, - only_obs = TRUE, - t_obs = c(2,4,6,8), - obs_type = c(1, 1, 1, 1), - seed = 123456, - res_var = ruv_multi, - t_init = 10 - ) - expect_equal(s_multi1$y, s_single$y) - - ## specifying ruv as multi-type when multiple obs_type - ## should now give different results - s_multi1 <- sim( - ode = pk1, - parameters = parameters, - n_ind = 1, - regimen = regimen, - only_obs = TRUE, - t_obs = c(2,4,6,8), - obs_type = c(1, 2, 1, 2), - seed = 123456, - res_var = ruv_multi, - t_init = 10 - ) - expect_equal(s_multi1$y[c(1, 3)], s_single$y[c(1,3)]) - expect_true(sum(abs(s_single$y[c(2,4)] - s_multi1$y[c(2,4)])) > 50) -}) - -test_that("multi-obs with baseline and obs_time = dose_time works correctly", { - tmp <- sim( - pk1, - parameters = parameters, - regimen = regimen, - t_obs = c(0, 0, 6, 6), - obs_type = c(1, 4, 1, 4), - only_obs = TRUE, - output_include = list(variables = TRUE), - return_design = F - ) - expect_equal( - tmp$y[tmp$obs_type == 1], - tmp$CONC[tmp$obs_type == 1] - ) - expect_equal( - tmp$y[tmp$obs_type == 4], - tmp$ACT[tmp$obs_type == 4] - ) -}) - -test_that("multi-obs model simulates correctly (bug 202205)", { - skip_on_cran() - - pars <- list( - CL = 23.9, - CLM = 5.19, - V = 107, - Q = 3.31, - V2 = 46.2, - VM = 111, - KA = 18.6, - TLAG = 0.415 - ) - covs <- list(WT = PKPDsim::new_covariate(80)) - mod <- new_ode_model( - code = "dAdt[0] = -KAi*A[0] \ - dAdt[1] = +KAi*A[0] - (CLi/Vi)*A[1] - (Qi/Vi)*A[1] + (Qi/V2i)*A[2] \ - dAdt[2] = +(Qi/Vi)*A[1] - (Qi/V2i)*A[2] \ - dAdt[3] = +(CLi/Vi)*A[1] - (CLMi/VMi)*A[3] \ - dAdt[4] = A[1]*1000.0/Vi \ - CONC = A[1]/(Vi/1000.0) \ - CONCM = A[3]/(VMi/1000.0) \ - CONCT = CONC+CONCM \ - ", - pk = "KAi = KA \ - CLi = CL * pow(WT/70.0, 0.75) \ - Vi = V *(WT/70.0) \ - Qi = Q * pow(WT/70.0, 0.75)\ - V2i = V2 * (WT/70.0) \ - CLMi = CLM * pow(WT/70.0, 0.75) \ - VMi = VM *(WT/70.0) \ - ", - parameters = pars, - declare_variables = c( "CLi", "Vi", "V2i", "Qi", "KAi", "VMi", "CLMi", "CONC", "CONCM", "CONCT" ), - obs = list( - variable = c( "CONC", "CONCM", "CONC", "CONCM" ) - ), - covariates = covs - ) - regimen <- new_regimen(amt = 3.7, n=10, interval = 24, type = "oral", cmt = 1) - data <- data.frame( - t = 70, - y = 1, - evid = 0, - loq = 0, - obs_type = 2, - dv = 1 - ) - t_obs <- seq(0, 50, .5) - pop <- sim( - ode = mod, - regimen = regimen, - parameters = pars, - covariates = covs, - t_obs = rep(t_obs, each = 4), - only_obs = TRUE, - obs_type = rep(1:4, length(t_obs)) - ) - # correct obs_types - expect_equal( - pop[pop$t %in% c(23, 47), ]$obs_type, - c(1,2,3,4,1,2,3,4) - ) - # correct output values - expect_equal( - round(pop[pop$t %in% c(23, 47), ]$y, 2), - c(0.52, 12.34, 0.52, 12.34, 0.63, 17.17, 0.63, 17.17) - ) - -}) diff --git a/tests/testthat/test_multiple_regimens.R b/tests/testthat/test_multiple_regimens.R deleted file mode 100644 index 59a312b1..00000000 --- a/tests/testthat/test_multiple_regimens.R +++ /dev/null @@ -1,24 +0,0 @@ -test_that("multiple regimens for multiple individuals can be simulated", { - # Uses a model defined in `setup.R` - cov_table <- data.frame(WT = rnorm(10, 70, 5)) - multi_regs <- list() - for(i in seq(cov_table$WT)) { - multi_regs[[i]] <- new_regimen(amt = 10 * cov_table$WT[i], interval = 12, type = "infusion") - } - class(multi_regs) <- "regimen_multiple" - - par <- list(CL = 5, V = 50) - reg <- new_regimen(amt = 2000, interval = 24, type = "infusion") - covariates = list(WT = new_covariate(1)) - res <- sim_ode( - ode = mod_1cmt_iv, - parameters = par, - covariates = covariates, - regimen = multi_regs, - only_obs = TRUE - ) - expect_equal(length(unique(res$id)),10) - expect_equal(sum(is.na(res$y)), 0) - # IDs ordered correctly: - expect_true(length(unique(res$id)) == 10 && all(diff(res$id) >= 0)) -}) diff --git a/tests/testthat/test_adherence.R b/tests/testthat/test_new_adherence.R similarity index 100% rename from tests/testthat/test_adherence.R rename to tests/testthat/test_new_adherence.R diff --git a/tests/testthat/test_new_ode_model.R b/tests/testthat/test_new_ode_model.R index 5b5141a6..0b3c101f 100644 --- a/tests/testthat/test_new_ode_model.R +++ b/tests/testthat/test_new_ode_model.R @@ -163,3 +163,108 @@ test_that("specifying overlapping covariates and variables throws error", { ) }) }) + +test_that("IOV: Change in F in 2nd bin is applied in 2nd bin and not later.", { + skip_on_cran() + # Previously this was an issue because F, when defined in pk_code(), was not updated before the + # dose was applied to the state vector, so the bioavailability was not applied at the right time. + # This was fixed by rearranging the order of execution in sim.cpp in the main loop. + + pars_iov_f <- list( + "CL" = 5, + "V" = 50, + "KA" = 1, + "F" = 1 + ) + pk_iov_f <- new_ode_model( + code = " + dAdt[1] = -KA * A[1] + dAdt[2] = +KA * A[1] -(CLi/V) * A[2] + ", + pk_code = " + Fi = F * exp(kappa_F); + CLi = CL; + ", + iov = list( + cv = list(F = 0.2), + n_bins = 3 + ), + obs = list(cmt = 2, scale = "V"), + dose = list(cmt = 1, bioav = "Fi"), + declare_variables = c("kappa_F", "Fi", "CLi"), + parameters = names(pars_iov_f), + cpp_show_code = FALSE + ) + reg <- new_regimen(amt = 800, interval = 24, n = 10, type = "oral") + + # For a first simulation, we're simulating with no variability across the IOV bins: + pars_iov_f$kappa_F_1 <- 0 + pars_iov_f$kappa_F_2 <- 0 + pars_iov_f$kappa_F_3 <- 0 + args_sim1 <- args <- list( + ode = pk_iov_f, + parameters = pars_iov_f, + regimen = reg, + only_obs = TRUE, + t_obs = seq(0, 50, .25), + iov_bins = c(0L, 24L, 48L, 9999L) + ) + # For a second simulation, we're applying a change in parameter for the 2nd bin (24-48 hrs). + # This should affect predictions from 24 onward. + pars_iov_f$kappa_F_2 <- 1 # 2nd bin + args_sim2 <- args <- list( + ode = pk_iov_f, + parameters = pars_iov_f, + regimen = reg, + only_obs = TRUE, + t_obs = seq(0, 50, .25), + iov_bins = c(0L, 24L, 48L, 9999L) + ) + res1 <- do.call("sim_ode", args = args_sim1) + res2 <- do.call("sim_ode", args = args_sim2) + expect_true(min(res1[res1$y != res2$y,]$t) <= 25) +}) + +describe("Models with bioavailability", { + parameters <- list(KA = 0.5, CL = 5, V = 50) + covs <- list(WT = new_covariate(50)) + reg <- new_regimen(amt = 100, n = 1, interval = 12, type = "bolus") + mod1 <- oral_1cmt_allometric # defined in setup.R using same covs and pars as above + y1 <- sim_ode(mod1, parameters = parameters, regimen = reg, covariates = covs, only_obs=TRUE)$y + + test_that("bioav numeric option working", { + skip_on_cran() + mod2 <- new_ode_model( + code = " + CLi = CL * pow(WT/70, 0.75) + dAdt[1] = -KA * A[1] + dAdt[2] = KA*A[1] - (CLi/V)*A[2] + ", + dose = list(cmt = 1, bioav = 0.5), + covariates = covs, + declare_variables = "CLi", + parameters = parameters + ) + + y2 <- sim_ode(mod2, parameters = parameters, regimen = reg, covariates = covs, only_obs=TRUE)$y + + expect_equal(round(y1,1), round(y2*2, 1)) + }) + + test_that("bioav string option working", { + skip_on_cran() + mod3 <- new_ode_model( + code = " + CLi = CL * pow(WT/70, 0.75) + dAdt[1] = -KA * A[1] + dAdt[2] = KA*A[1] - (CLi/V)*A[2] + ", + dose = list(cmt = 1, bioav = "WT/70"), + covariates = covs, + declare_variables = "CLi", + parameters = parameters + ) + y3 <- sim_ode(mod3, parameters = parameters, regimen = reg, covariates = covs, only_obs=TRUE)$y + expect_equal(round(y1*(50/70),1), round(y3, 1)) + }) +}) \ No newline at end of file diff --git a/tests/testthat/test_parameters_table.R b/tests/testthat/test_parameters_table.R index 1f9af030..4dcacdd1 100644 --- a/tests/testthat/test_parameters_table.R +++ b/tests/testthat/test_parameters_table.R @@ -1,23 +1,16 @@ test_that("Simulating with table of parameters works", { skip_on_cran() parameters_table <- data.frame( + KA = rnorm(10, 1, .1), CL = rnorm(10, 5, 5), V = rnorm(10, 5, 0.5) ) - mod <- new_ode_model( - code = " - dAdt[1] = -(CL/V) * A[1]; - ", - obs = list(cmt = 1, scale = "V"), - covariates = list(WT = new_covariate(70)), - cpp_show_code=FALSE - ) par <- list(CL = 5, V = 50) reg <- new_regimen(amt = 2000, interval = 24, type = "infusion") covariates = list(WT = new_covariate(1)) res <- sim_ode( - ode = mod, + ode = oral_1cmt_allometric, parameters_table = parameters_table, covariates = covariates, regimen = reg, diff --git a/tests/testthat/test_nlmixr_translation.R b/tests/testthat/test_pkpdsim_to_nlmixr.R similarity index 100% rename from tests/testthat/test_nlmixr_translation.R rename to tests/testthat/test_pkpdsim_to_nlmixr.R diff --git a/tests/testthat/test_reparametrization.R b/tests/testthat/test_reparametrization.R index 747d3734..f3352def 100644 --- a/tests/testthat/test_reparametrization.R +++ b/tests/testthat/test_reparametrization.R @@ -1,3 +1,9 @@ +# Uses model and covariates defined in setup.R (conditional, NOT_CRAN only): +# - model_carreno +# - covs_carreno + +skip_on_cran() + par_orig <- list( V = 25.76, SCLSlope = 0.036, @@ -6,33 +12,8 @@ par_orig <- list( SCLInter = 0.18, TDM_INIT = 0) -covs <- list( - CRCL = PKPDsim::new_covariate(5), - CL_HEMO = PKPDsim::new_covariate(0) -) -model <- new_ode_model( # Carreno et al - code = " - CLi = SCLInter + SCLSlope * (CRCL*16.667) + CL_HEMO \ - Vi = V \ - Qi = K12 * Vi \ - V2i = Qi / K21 \ - dAdt[0] = -(CLi/V)*A[0] - K12*A[0] + K21*A[1] \ - dAdt[1] = K12*A[0] - K21*A[1] \ - dAdt[2] = A[0]/V - ", - parameters = par_orig, - declare_variables = c("CLi", "Qi", "Vi", "V2i"), - covariates= covs, - obs = list(cmt = 1, scale = "V"), - reparametrization = list( - "CL" = "SCLInter + SCLSlope * (CRCL*16.667) + CL_HEMO", - "V" = "V", - "Q" = "K12 * V", - "V2" = "(K12 * V) / K21" - ) -) -pars_covs_comb <- join_cov_and_par(covs, par_orig) -pars <- reparametrize(model, pars_covs_comb, covariates = covs) +pars_covs_comb <- join_cov_and_par(covs_carreno, par_orig) +pars <- reparametrize(model_carreno, pars_covs_comb, covariates = covs_carreno) reg <- new_regimen( amt = 1000, @@ -44,10 +25,10 @@ reg <- new_regimen( n_ss = 50 ) s <- sim( - ode = model, + ode = model_carreno, parameters = par_orig, regimen = reg, - covariates = covs, + covariates = covs_carreno, t_obs = c(0,24) ) @@ -96,4 +77,3 @@ test_that("Reparametrized model and analytics equations match", { cmin_ss_ode <- s[s$comp == "obs",]$y[1] expect_equal(cmin_ss_ode, cmin_ss_lin) }) - diff --git a/tests/testthat/test_sim.R b/tests/testthat/test_sim.R index cb676754..62a29a1e 100644 --- a/tests/testthat/test_sim.R +++ b/tests/testthat/test_sim.R @@ -578,3 +578,995 @@ test_that("t_max is shifted correctly when t_ss != 0", { )) expect_equal(sum(is.na(evtab3)), 0) }) + +describe("IOV", { + + reg_iov <- new_regimen( + amt = 100, + interval = 24, + n = 5, + type = "infusion" + ) + iov_var <- 0.3 ^ 2 # 30% IOV + + test_that("Throws error when `iov_bins` supplied but not present in model", { + expect_error({ + sim( + ode = pk_iov_none, + parameters = pars_iov_no_iov, + regimen = reg_iov, + omega = c( + 0.3 # IIV in CL + ), + n = 1, + iov_bins = c(0, 24, 48, 72, 999), # !! + omega_type = "normal", + only_obs = TRUE, + output_include = list(parameters = TRUE, variables = TRUE) + ) + }, "No IOV implemented for this model") + }) + + test_that("Throws error when number of `iov_bins` is higher than allowed for model", { + expect_error({ + sim( + ode = pk_iov, + parameters = pars_iov, + regimen = reg_iov, + omega = c( + iov_var, # IOV in CL + 0, iov_var, + 0, 0, iov_var, + 0, 0, 0, iov_var, + 0, 0, 0, 0, 0.3 # IIV in CL + ), + n = 1, + iov_bins = c(0, 24, 48, 72, 999), # one bin too many + omega_type = "normal", + only_obs = TRUE, + output_include = list(parameters = TRUE, variables = TRUE) + ) + }, "Number of allowed IOV bins for model is lower") + }) + + test_that("Throws warning when number of `iov_bins` is lower than allowed for model", { + expect_warning({ + sim( + ode = pk_iov, + parameters = pars_iov, + regimen = reg_iov, + omega = c( + iov_var, # IOV in CL + 0, iov_var, + 0, 0, iov_var, + 0, 0, 0, iov_var, + 0, 0, 0, 0, 0.3 # IIV in CL + ), + n = 1, + iov_bins = c(0, 24, 999), # one bin too few + omega_type = "normal", + only_obs = TRUE, + output_include = list(parameters = TRUE, variables = TRUE) + )}, "Number of allowed IOV bins for model is higher" + ) + }) + + test_that("IOV is added to parameters", { + skip_on_cran() + set.seed(32) + + dat <- sim( + ode = pk_iov, + parameters = pars_iov, + regimen = reg_iov, + omega = c( + iov_var, # IOV in CL + 0, iov_var, + 0, 0, iov_var, + 0, 0, 0, iov_var, + 0, 0, 0, 0, 0.3 # IIV in CL + ), + n = 5, + iov_bins = c(0, 24, 48, 72), + omega_type = "normal", + only_obs = TRUE, + output_include = list(parameters = TRUE, variables = TRUE) + ) + + expect_equal( + signif(dat$kappa_CL[dat$t == 12], 4), + signif(dat$kappa_CL_1[dat$t == 1], 4) + ) + expect_equal( + signif(dat$kappa_CL[dat$t == 36], 4), + signif(dat$kappa_CL_2[dat$t == 1], 4) + ) + expect_equal( + signif(dat$kappa_CL[dat$t == 48], 4), + signif(dat$kappa_CL_3[dat$t == 1], 4) + ) + expect_equal( + signif(exp(dat$kappa_CL + dat$eta_CL) * dat$CL, 4), + signif(dat$CL_iov, 4) + ) + }) + + test_that("error is not invoked when using parameters_table", { + parameters_table <- data.frame( + CL = runif(10), V = runif(10), KA = runif(10), eta_CL = runif(10), + kappa_CL_1 = 0, kappa_CL_2 = 0, kappa_CL_3 = 0, kappa_CL_4 = 0 + ) + + # specifying both parameters_table but for a model with IOV should not fail! + expect_silent( + dat <- sim( + ode = pk_iov, + parameters_table = parameters_table, + regimen = reg_iov, + omega = c( + iov_var, # IOV in CL + 0, iov_var, + 0, 0, iov_var, + 0, 0, 0, iov_var, + 0, 0, 0, 0, 0.3 # IIV in CL + ), + n = 5, + iov_bins = c(0, 24, 48, 72), + omega_type = "normal", + only_obs = TRUE, + output_include = list(parameters = TRUE, variables = TRUE) + ) + ) + + # specifying both parameters and parameters_table should fail + expect_error( + dat <- sim( + ode = pk_iov, + parameters = pars_iov, + parameters_table = parameters_table, + regimen = reg_iov, + omega = c( + iov_var, # IOV in CL + 0, iov_var, + 0, 0, iov_var, + 0, 0, 0, iov_var, + 0, 0, 0, 0, 0.3 # IIV in CL + ), + n = 5, + iov_bins = c(0, 24, 48, 72), + omega_type = "normal", + only_obs = TRUE, + output_include = list(parameters = TRUE, variables = TRUE) + ) + ) + }) +}) + +describe("Compare results from sims with references", { + + # Uses models defined in setup.R: + # - mod_1cmt_oral + # - mod_1cmt_iv + # - dose_in_cmt_2 + + pk1cmt_oral_anal = function(t, dose, KA, V, CL) { + dose*KA/(V*(KA-CL/V))*(exp(-(CL/V) * t)-exp(-KA * t)) + } + + pk1cmt_oral_code <- new_ode_model( + code = "dAdt[1] = -KA*A[1]; dAdt[2] = KA*A[1] - (CL/V)*A[2];", + obs = list(cmt = 2, scale = "V") + ) + + test_that("Model from library and custom C++ and code matches analytic solution", { + p <- list(KA = 1, CL = 5, V = 50) + t_obs <- c(0:72) + t_obs2 <- t_obs + 0.1234 # also needs to be producing results with non-integer times + dose <- 100 + t_dose <- c(0) + regimen <- new_regimen(amt=dose, times = t_dose, type = "oral") + + pk1cmt_oral_lib <- sim_ode( + ode = mod_1cmt_oral, + parameters = p, + regimen = regimen, + t_obs = t_obs, + int_step_size = 0.1, + duplicate_t_obs = TRUE, + only_obs=TRUE + ) + + pk1cmt_oral_code_res <- sim_ode( + ode = pk1cmt_oral_code, + parameters = p, + duplicate_t_obs = TRUE, + regimen=regimen, + t_obs=t_obs, + int_step_size = 0.1, + only_obs=TRUE + ) + + pk1cmt_oral_anal_res <- pk1cmt_oral_anal(t_obs, dose, p$KA, p$V, p$CL) + expect_equal(round(pk1cmt_oral_lib$y, 3), round(pk1cmt_oral_anal_res, 3)) + expect_equal(round(pk1cmt_oral_code_res$y, 3), round(pk1cmt_oral_anal_res, 3)) + }) + + + test_that("precision in time does not impact # obs returned", { + regimen_mult <- new_regimen( + amt = rep(12.8, 3), + times = c(0, 6, 12), + type="infusion", + t_inf = 2 + ) + t_obs <- c(11.916, 14.000, 16.000, 17.000, 30) + tmp <- sim_ode( + ode = mod_1cmt_iv, + parameters = list(CL = 5, V = 50), + regimen = regimen_mult, + t_obs = t_obs, + only_obs = TRUE + ) + expect_equal(tmp$t, t_obs) + }) + + test_that("test bug EmCo 20150925", { + xtim <- c(0, 2, 4, 8, 12, 24) + sujdos <- 320 + param <- list(KA = 1.8, V = 30, CL = 1.7) + regim <- new_regimen(amt = sujdos, times = c(0, 12), type= "bolus") + out <- sim_ode(ode = mod_1cmt_oral, parameters=param, regimen=regim, t_obs = xtim, only_obs = TRUE) + expect_equal(out$t, xtim) + }) + + test_that("model size is appropriate (bug: JeHi 20151204)", { + skip_on_cran() + pk3cmt <- new_ode_model( + code = " + dAdt[1] = -KA*A[1]; + dAdt[2] = KA*A[1] -(Q/V)*A[2] + (Q/V2)*A[3] -(CL/V)*A[2]; + dAdt[3] = -(Q/V2)*A[3] + (Q/V)*A[2]; + ", + obs = list(cmt = 2, scale = "V") + ) + expect_equal( attr(pk3cmt, "size"), 3) + }) + + test_that("Dose is added to correct compartment: specified by model", { + set.seed(90) + p <- list(CL = 1, V = 10, KA = 0.5, S2 = .1) + r <- new_regimen(amt = 100, times = c(0), type = "infusion") + dat <- sim_ode( + ode = dose_in_cmt_2, + n_ind = 1, + omega = cv_to_omega(par_cv = list("CL"=0.1, "V"=0.1, "KA" = .1), p), + parameters = p, + regimen = r, + verbose = FALSE, + t_max = 48 + ) + # Dose should be in cmt 2 + expect_equal(dat$y[dat$comp == 1], rep(0, 50)) + expect_true(all(dat$y[dat$comp == 2][-1] > 0)) + }) + + test_that("Dose is added to correct compartment: override model by regimen", { + set.seed(60) + p <- list(CL = 1, V = 10, KA = 0.5, S2 = .1) + r <- new_regimen( + amt = c(100, 100, 100), + times = c(0, 6, 12), + cmt = c(1,2,3), + type = "bolus" + ) + dat <- sim_ode( + ode = dose_in_cmt_2, + n_ind = 1, + omega = cv_to_omega(par_cv = list("CL"=0.1, "V"=0.1, "KA" = .1), p), + parameters = p, + regimen = r, + verbose = FALSE, + t_max = 48 + ) + # Dose should be in cmt 1, 2 and 3 + expect_true(all(dat$y[dat$comp == 1 & dat$t > 0] > 0)) + expect_true(max(diff(dat$y[dat$comp == 2])) > 95) + expect_true(max(diff(dat$y[dat$comp == 3])) > 95) + }) + + test_that("Infusion works for all compartments", { + set.seed(44) + # Part 1: Specify cmt only with model + p <- list(CL = 1, V = 10, KA = 0.5, S2 = .1) + r <- new_regimen( + amt = c(100, 100, 100), + times = c(0, 6, 12), + cmt = c(1,2,3), + t_inf = 3, + type = "infusion" + ) + dat <- sim_ode( + ode = dose_in_cmt_2, + n_ind = 1, + omega = cv_to_omega(par_cv = list("CL"=0.1, "V"=0.1, "KA" = .1), p), + parameters = p, + regimen = r, + verbose = FALSE, + t_max = 48 + ) + expect_true(all(dat$y[dat$comp == 1 & dat$t > 0 ] > 0)) + expect_true(max(diff(dat$y[dat$comp == 2])) > 25) + expect_true(max(diff(dat$y[dat$comp == 3])) > 25) + expect_equal(round(max(dat$y[dat$comp == 2]), 1), 131.2) + expect_equal(round(max(dat$y[dat$comp == 3]), 1), 148.4) + }) + + test_that("Duplicate obs returned when specified in arg", { + # for first 2 doses, infusion time will just be ignored, but a value has to be specified in the vector + p <- list(CL = 1, V = 10, KA = 0.5, S2=.1) + r <- new_regimen( + amt = c(100, 100, 100, 100), + times = c(0, 6, 12, 18), + cmt = c(2, 2, 1, 1), + t_inf = c(1, 1, 1, 1), + type = c("bolus", "bolus", "infusion", "infusion") + ) + dat <- sim_ode( + ode = mod_1cmt_oral, + n_ind = 1, + omega = cv_to_omega(par_cv = list("CL"=0.1, "V"=0.1, "KA" = .1), p), + parameters = p, + regimen = r, + t_obs = c(1, 2, 3, 4, 4, 4, 6), ## see duplicate obs here + duplicate_t_obs = T, + only_obs = FALSE + ) + expect_equal(length(dat[dat$t == 4,]$y), 9) + expect_equal(length(dat$y), 21) + expect_equal(sum(is.na(dat$y)), 0) + }) + + test_that("Custom t_obs is returned", { + t_obs <- seq(from = 0, to = 24, by = .1) + p <- list(CL = 1, V = 10, KA = 0.5, S2=.1) + r <- new_regimen( + amt = c(100, 100, 100, 100), + times = c(0, 6, 12, 18), + cmt = c(2, 2, 1, 1), + t_inf = c(1, 1, 1, 1), + type = c("bolus", "bolus", "infusion", "infusion") + ) + dat <- sim_ode( + ode = mod_1cmt_oral, + n_ind = 1, + omega = cv_to_omega(par_cv = list("CL"=0.1, "V"=0.1, "KA" = .1), p), + parameters = p, + regimen = r, + t_obs = t_obs, + only_obs = T + ) + expect_equal(mean(diff(t_obs)), mean(diff(dat$t))) + }) + + test_that("if covariate time is at end of infusion, end of infusion is still recorded", { + # Bug reported by JF + pop_est <- list(CL = 1.08, V = 0.98) + regimen <- new_regimen( + amt = c(1500, 1000, 1500, 1500, 1500, 1500, 1500), + type = "infusion", + t_inf = c(2, 1, 2, 2, 1, 1, 1), + times = c(0, 10.8666666666667, 20.4333333333333, 32.0666666666667, 46.9, 54.9, 62.9 ) + ) + covs <- list( + WT = new_covariate(value = c(60, 65), times = c(0, 47.9)), + CRCL = new_covariate(8), CVVH = new_covariate(0) + ) + pksim <- sim( + ode = mod_1cmt_iv, + parameters = pop_est, + covariates = covs, + regimen = regimen, + checks = TRUE, + only_obs = TRUE + ) + expect_true(all(pksim$y < 1000)) + }) + + test_that("Covariate table simulation runs", { + # this test used to be in the covariate_table_to_list file but + # makes more sense here. + p <- list(CL = 5, V = 50) + reg <- new_regimen (amt = 100, n = 4, interval = 12, type = "bolus", cmt=1) + om <- c(0.01, 1, 0.01) + cov_table <- data.frame( + id=c(1, 1, 2, 3), + WT = c(40, 45, 50, 60), + SCR = c(50, 150, 90,110), + t = c(0, 5, 0, 0) + ) + + dat <- sim( + mod_1cmt_iv, + parameters = p, + regimen = reg, + covariates_table = cov_table, + covariates_implementation = list(SCR = "interpolate"), + omega = NULL, + n_ind = 3, + only_obs = T, + output_include = list(parameters = TRUE, covariates=TRUE) + ) + expect_equal(length(unique(dat$id)), 3) + }) + +}) + +describe("Simulate with multiple observation types", { + + # Uses model defined in setup.R: + # - pk_multi_obs + # - vars_multi_obs + + regimen_multi_obs <- new_regimen(amt = 100, interval = 12, n = 5, type="infusion", t_inf = 1) + parameters_multi_obs <- list("CL" = 15, "V" = 150) + omega_multi_obs <- PKPDsim::cv_to_omega(list("CL" = 0.2, "V" = 0.2), parameters_multi_obs[1:2]) + + test_that("obs types are output by `sim`", { + obs_type <- c(1,2,1,3,1) + data <- sim( + ode = pk_multi_obs, + parameters = list(CL = 20, V = 200), + regimen = regimen_multi_obs, + int_step_size = 0.1, + only_obs = TRUE, + obs_type = obs_type, + t_obs = c(2, 4, 6, 8, 12), + output_include = list("variables" = TRUE) + ) + expect_equal(data$obs_type, obs_type) + expect_equal(data$y, diag(as.matrix(data[1:5,5+obs_type]))) + }) + + + test_that("check obs at same timepoint but with different obs_type", { + obs_type <- c(1,2,1,3,1) + t_same <- sim( + ode = pk_multi_obs, + parameters = list(CL = 20, V = 200), + regimen = regimen_multi_obs, + int_step_size = 0.1, + only_obs = TRUE, + obs_type = obs_type, + t_obs = c(2, 4, 4, 8, 8), + output_include = list("variables" = TRUE) + ) + expect_equal(t_same$t[2], t_same$t[3]) + expect_equal(t_same$t[4], t_same$t[5]) + expect_equal(t_same$y[3], t_same$METAB[3]) + expect_equal(t_same$y[2], t_same$CONC[2]) + expect_equal(t_same$y[5], t_same$METAB2[5]) + expect_equal(t_same$y[4], t_same$CONC[4]) + }) + + + test_that("check that residual error correctly applied to right var", { + obs_type <- c(1,2,1,3,1) + set.seed(12345) + ruv_term3 <- list(prop = c(0, 0, 0.1), add = c(0, 0, 0.1)) + ruv_term1 <- list(prop = c(.1, 0, 0), add = c(1, 0, 0)) + + data2 <- sim( + ode = pk_multi_obs, + parameters = list(CL = 20, V = 200), + regimen = regimen_multi_obs, + int_step_size = 0.1, + only_obs = TRUE, + obs_type = obs_type, + t_obs = c(2, 4, 6, 8, 12), + output_include = list("variables" = TRUE), + res_var = ruv_term3 + ) + y <- diag(as.matrix(data2[1:5, 5 + obs_type])) + + # no residual error for obs_type 1 and 2, only 3 + expect_equal(data2$y[-4], y[-4]) + expect_true(data2$y[-4][4] != y[4]) + + data3 <- sim( + ode = pk_multi_obs, + parameters = list(CL = 20, V = 200), + regimen = regimen_multi_obs, + int_step_size = 0.1, + only_obs = TRUE, + obs_type = obs_type, + t_obs = c(2, 4, 6, 8, 12), + output_include = list("variables" = TRUE), + res_var = ruv_term1 + ) + y <- diag(as.matrix(data2[1:5, 5 + obs_type])) + + # no residual error for obs_type 2/3, only and 1 + expect_equal(data3$y[-c(1,3,5)], y[-c(1,3,5)]) + expect_true(all(data3$y[c(1,3,5)] != y[c(1,3,5)])) + }) + + + test_that("specifying ruv as multi-type when only 1 obs_type", { + obs_type <- c(1,2,1,3,1) + ruv_single <- list(prop = 0.1, add = 1) + ruv_multi <- list(prop = c(0.1, 1), add = c(1, 20)) + + s_single <- sim( + ode = pk_multi_obs, + parameters = parameters_multi_obs, + n_ind = 1, + regimen = regimen_multi_obs, + only_obs = TRUE, + t_obs = c(2,4,6,8), + seed = 123456, + res_var = ruv_single + ) + + ## specified as multi, but obs_type is all 1, so should + ## give same results as s_single + s_multi1 <- sim( + ode = pk_multi_obs, + parameters = parameters_multi_obs, + n_ind = 1, + regimen = regimen_multi_obs, + only_obs = TRUE, + t_obs = c(2,4,6,8), + obs_type = c(1, 1, 1, 1), + seed = 123456, + res_var = ruv_multi, + t_init = 10 + ) + expect_equal(s_multi1$y, s_single$y) + + ## specifying ruv as multi-type when multiple obs_type + ## should now give different results + s_multi1 <- sim( + ode = pk_multi_obs, + parameters = parameters_multi_obs, + n_ind = 1, + regimen = regimen_multi_obs, + only_obs = TRUE, + t_obs = c(2,4,6,8), + obs_type = c(1, 2, 1, 2), + seed = 123456, + res_var = ruv_multi, + t_init = 10 + ) + expect_equal(s_multi1$y[c(1, 3)], s_single$y[c(1,3)]) + expect_true(sum(abs(s_single$y[c(2,4)] - s_multi1$y[c(2,4)])) > 50) + }) + + test_that("multi-obs with baseline and obs_time = dose_time works correctly", { + tmp <- sim( + pk_multi_obs, + parameters = parameters_multi_obs, + regimen = regimen_multi_obs, + t_obs = c(0, 0, 6, 6), + obs_type = c(1, 4, 1, 4), + only_obs = TRUE, + output_include = list(variables = TRUE), + return_design = F + ) + expect_equal( + tmp$y[tmp$obs_type == 1], + tmp$CONC[tmp$obs_type == 1] + ) + expect_equal( + tmp$y[tmp$obs_type == 4], + tmp$ACT[tmp$obs_type == 4] + ) + }) + + test_that("multi-obs model simulates correctly (bug 202205)", { + skip_on_cran() + + pars <- list( + CL = 23.9, + CLM = 5.19, + V = 107, + Q = 3.31, + V2 = 46.2, + VM = 111, + KA = 18.6, + TLAG = 0.415 + ) + covs <- list(WT = PKPDsim::new_covariate(80)) + mod <- new_ode_model( + code = "dAdt[0] = -KAi*A[0] \ + dAdt[1] = +KAi*A[0] - (CLi/Vi)*A[1] - (Qi/Vi)*A[1] + (Qi/V2i)*A[2] \ + dAdt[2] = +(Qi/Vi)*A[1] - (Qi/V2i)*A[2] \ + dAdt[3] = +(CLi/Vi)*A[1] - (CLMi/VMi)*A[3] \ + dAdt[4] = A[1]*1000.0/Vi \ + CONC = A[1]/(Vi/1000.0) \ + CONCM = A[3]/(VMi/1000.0) \ + CONCT = CONC+CONCM \ + ", + pk = "KAi = KA \ + CLi = CL * pow(WT/70.0, 0.75) \ + Vi = V *(WT/70.0) \ + Qi = Q * pow(WT/70.0, 0.75)\ + V2i = V2 * (WT/70.0) \ + CLMi = CLM * pow(WT/70.0, 0.75) \ + VMi = VM *(WT/70.0) \ + ", + parameters = pars, + declare_variables = c( "CLi", "Vi", "V2i", "Qi", "KAi", "VMi", "CLMi", "CONC", "CONCM", "CONCT" ), + obs = list( + variable = c( "CONC", "CONCM", "CONC", "CONCM" ) + ), + covariates = covs + ) + regimen <- new_regimen(amt = 3.7, n=10, interval = 24, type = "oral", cmt = 1) + data <- data.frame( + t = 70, + y = 1, + evid = 0, + loq = 0, + obs_type = 2, + dv = 1 + ) + t_obs <- seq(0, 50, .5) + pop <- sim( + ode = mod, + regimen = regimen, + parameters = pars, + covariates = covs, + t_obs = rep(t_obs, each = 4), + only_obs = TRUE, + obs_type = rep(1:4, length(t_obs)) + ) + # correct obs_types + expect_equal( + pop[pop$t %in% c(23, 47), ]$obs_type, + c(1,2,3,4,1,2,3,4) + ) + # correct output values + expect_equal( + round(pop[pop$t %in% c(23, 47), ]$y, 2), + c(0.52, 12.34, 0.52, 12.34, 0.63, 17.17, 0.63, 17.17) + ) + + }) + +}) + +describe("Multiple regimens", { + test_that("multiple regimens for multiple individuals can be simulated", { + + # Uses a model defined in `setup.R` + cov_table <- data.frame(WT = rnorm(10, 70, 5)) + multi_regs <- list() + for(i in seq(cov_table$WT)) { + multi_regs[[i]] <- new_regimen(amt = 10 * cov_table$WT[i], interval = 12, type = "infusion") + } + class(multi_regs) <- "regimen_multiple" + + par <- list(CL = 5, V = 50) + reg <- new_regimen(amt = 2000, interval = 24, type = "infusion") + covariates = list(WT = new_covariate(1)) + res <- sim_ode( + ode = mod_1cmt_iv, + parameters = par, + covariates = covariates, + regimen = multi_regs, + only_obs = TRUE + ) + expect_equal(length(unique(res$id)),10) + expect_equal(sum(is.na(res$y)), 0) + # IDs ordered correctly: + expect_true(length(unique(res$id)) == 10 && all(diff(res$id) >= 0)) + }) + +}) + + +describe("Mixture models", { + # Uses model and covariates defined in setup.R: + # - mod_mixture + # - covs_mixture + + par_mixture <- list(CL = 3, V = 50) + reg_mixture <- new_regimen(amt = 250, n = 5, interval = 6, type = 'infusion', t_inf = 1) + t_obs_mixture <- seq(0, 36, 4) + + test_that("mixture model works properly for single patient", { + res0 <- sim_ode(mod_mixture, parameters = par_mixture, regimen = reg_mixture, covariates = covs_mixture, t_obs = t_obs_mixture, only_obs=T) # mixture_group not supplied + res1 <- sim(mod_mixture, parameters = par_mixture, regimen = reg_mixture, t_obs = t_obs_mixture, covariates = covs_mixture, mixture_group = 1, only_obs=T) + res2 <- sim(mod_mixture, parameters = par_mixture, regimen = reg_mixture, t_obs = t_obs_mixture, covariates = covs_mixture, mixture_group = 2, only_obs=T) + expect_equal(round(res0[res0$t == 24,]$y, 2), 9.07) # should use whatever is in `parameters` + expect_equal(round(res1[res1$t == 24,]$y, 2), 5.82) + expect_equal(round(res2[res2$t == 24,]$y, 2), 1.15) + }) + + test_that("mixture model works properly when vectorized (using parameters_table)", { + partab <- data.frame(CL = rep(0, 6), V = rep(50, 6)) + suppressMessages({ + expect_error(sim_ode(mod_mixture, parameters_table = partab, regimen = reg_mixture, t_obs = t_obs_mixture, covariates = covs_mixture, mixture_group = 1, only_obs=T)) + res1 <- sim(mod_mixture, parameters_table = partab, regimen = reg_mixture, t_obs = t_obs_mixture, covariates = covs_mixture, mixture_group = rep(1, 6), only_obs=T) + res2 <- sim(mod_mixture, parameters_table = partab, regimen = reg_mixture, t_obs = t_obs_mixture, covariates = covs_mixture, mixture_group = rep(c(1,2), 3), only_obs=T, output_include = list(parameters = TRUE)) + }) + expect_equal(round(res1[res1$t == 24,]$y, 2), rep(5.82, 6)) + expect_equal(round(res2[res2$t == 24,]$y, 2), rep(c(5.82, 1.15), 3)) + expect_equal(res2[res2$id == 1,]$CL[1], 5) + expect_equal(res2[res2$id == 2,]$CL[1], 15) + expect_equal(res2[res2$id == 3,]$CL[1], 5) + }) + + test_that("mixture model works properly when vectorized (using covariates_table)", { + covtab <- data.frame(ID = 1:8, WT = rep(seq(40, 130, 30), 2)) + suppressMessages({ + expect_error(sim(mod_mixture, parameters = par_mixture, covariates_table = covtab, regimen = reg_mixture, t_obs = t_obs_mixture, mixture_group = 1, only_obs=T)) + res <- sim(mod_mixture, parameters = par_mixture, covariates_table = covtab, regimen = reg_mixture, t_obs = t_obs_mixture, mixture_group = rep(c(1, 2), each=4), only_obs=T) + }) + expect_equal(round(res[res$t == 24,]$y, 2), c(9.39, 5.82, 3.83, 2.65, 2.99, 1.15, 0.52, 0.25)) + }) + +}) + +describe("Compartment mapping", { + + pk1cmt_oral_cmt_mapping <- new_ode_model( + code = "dAdt[1] = -KA*A[1]; dAdt[2] = KA*A[1] - (CL/V)*A[2];", + obs = list(cmt = 2, scale = "V"), + cmt_mapping = list(oral = 1, infusion = 2, bolus = 2) + ) + + test_that("Compartment mapping is added to attributes", { + expect_equal(attr(pk1cmt_oral_cmt_mapping, "cmt_mapping")[["oral"]], 1) + expect_equal(attr(pk1cmt_oral_cmt_mapping, "cmt_mapping")[["infusion"]], 2) + }) + + test_that("Admin route is interpreted and simulated correctly", { + regimen <- new_regimen( + amt = c(100, 100, 100, 100), + times = c(0, 12, 24, 36), + type = c("oral", "oral", "infusion", "infusion"), + t_inf = c(0, 0, 1, 1) + ) + p <- list(KA = 1, CL = 5, V = 50) + res <- sim_ode( + ode = pk1cmt_oral_cmt_mapping, + parameters = p, + regimen = regimen, + only_obs = FALSE + ) + expect_equal(round(res$y[res$comp == 1 & res$t == 25], 4), 2e-04) + expect_true(res$y[res$comp == 2 & res$t == 25] >= 100) + }) + + test_that("multiple scaling types on one compartment works", { + skip_on_cran() + mod <- new_ode_model( + code = " + dAdt[1] = -KA * A[1]; + dAdt[2] = -(CL/V) * A[2] + KA*A[1]; + ", + obs = list( + cmt = c(2, 2), + scale = c(1, "V"), + label = c("abs", "conc") + ), + cpp_show_code = FALSE + ) + par <- list(CL = 5, V = 50, KA = .5) + reg <- new_regimen(amt = 100, n = 5, interval = 12) + res <- sim_ode( + ode = mod, + parameters = par, + regimen = reg, + only_obs = TRUE + ) + dat <- cbind(res[res$comp == "abs",]$y, res[res$comp == "conc",]$y) + expect_true("PKPDsim" %in% class(mod)) + expect_equal(length(unique(res$comp)), 2) + expect_equal(round(dat[,1],1), round(dat[,2]*par$V,1)) + }) + +}) + +describe("Model simulation with lagtime", { + + test_that("dose dump after lagtime in correct order in output data", { + skip_on_cran() + reg <- new_regimen(amt = 500, n = 4, interval = 12, type = 'oral') + pars <- list(CL = 5, V = 50, KA = 0.5, TLAG = 0.83) + dat <- sim_ode( + ode = mod_1cmt_oral_lagtime, + regimen = reg, + parameters = pars, + only_obs = FALSE + ) + ## Change after RXR-2394: time of TLAG not in dataset anymore unless requested by user in t_obs + ## Before change: expect_equal(round(dat[dat$t == 12.83 & dat$comp == 1,]$y, 1), c(1.2, 501.2)) + ## After change: + expect_equal(nrow(dat[dat$t == 12.83,]), 0) + ## When grid requested by user, lagtime should be visible + dat <- sim_ode( + ode = mod_1cmt_oral_lagtime, + regimen = reg, + parameters = pars, + t_obs = seq(0, 1, 0.01), + only_obs = TRUE + ) + tmp <- dat[dat$t >= 0.82 & dat$t <= 0.85, ] + expect_equal(tmp$t, c(0.82, 0.83, 0.84, 0.85)) + expect_equal(round(tmp$y, 2), c(0, 0, 0.05, 0.10)) + }) + +}) + + +describe("Models with t_init", { + # Uses model defined in setup.R (conditional, NOT_CRAN only): + # - mod_t_init + + skip_on_cran() + + # Test t_init functionality + ## e.g. TDM before first dose: + ## at t=-8, conc=10000 + ## Use this as true value in the simulations + par_t_init <- list(CL = 7.67, V = 97.7, TDM_INIT = 500) + reg_t_init <- new_regimen(amt = 100000, times=c(0, 24), type="bolus") + s <- sim( + ode = mod_t_init, + parameters = par_t_init, + n_ind = 1, + regimen = reg_t_init, + only_obs = TRUE, + output_include = list(variables = TRUE), + t_init = 10 + ) + + test_that("TDM before first dose is considered a true initial value", { + expect_equal(sum(is.na(s$y)), 0) + expect_equal(s$y[1], 500) + expect_equal(round(s$y[3]), 427) + expect_equal(round(s$y[13]),1157) + }) + + test_that("Variables are set (also in first row) when TDM before first dose", { + expect_equal(round(s$CONC[1:5], 1), c(500, 462.2, 427.3, 395.1, 365.3)) + }) + +}) + +describe("Time rounding issues", { + + # rounding time should not produce NAs in sim + ## time-rounding bug 20170804 + + test_that("No NAs related to rounding", { + # Uses models defined in setup.R + p <- list(CL = 5, V = 50, Q = 10, V2 = 150) + r1 <- new_regimen(amt = 100, times = c(0, 24, 36), type = "infusion") + dat1 <- sim_ode( + ode = mod_2cmt_iv, + parameters = p, + regimen = r1, + t_obs=seq(0, 150, length.out = 100) + ) + expect_equal(sum(is.na(dat1$y)), 0) + }) + +}) + +describe("Timevarying covariates are handled properly", { + # Uses model defined in setup.R (conditional, NOT_CRAN only): + # - mod_2cmt_timevar + + skip_on_cran() + + par_timevar <- list(CL = 3, V = 50, Q = 2.5, V2 = 70) + reg_timevar <- new_regimen( + amt = 250, + n = 60, + interval = 6, + type = 'infusion', + t_inf = 1 + ) + + test_that("timevarying covariates handled", { + # CLi changes by several orders of magnitude after + # steady state is achieved, which should produce + # a new steady state that is much lower + covs <- list( + CRCL = new_covariate( + value = c(4.662744, 5.798767, 6.195943, 6.2, 600), + times = c(0, 18, 29, 209, 210), + implementation = "locf" + ) + ) + t_obs <- seq(0, 360, 0.1) + sim1 <- sim_ode( + mod_2cmt_timevar, + parameters = par_timevar, + regimen = reg_timevar, + covariates = covs, + only_obs = TRUE, + t_obs = t_obs, + output_include = list(parameters = TRUE, covariates = TRUE) + ) + expect_equal(sim1$CRCL[sim1$t == 209], 6.2) + expect_equal(sim1$CRCL[sim1$t == 210], 600) + expect_true(sim1$y[sim1$t == 35 * 6] > 10 * sim1$y[sim1$t == 60 * 6]) + }) + + test_that("timevarying covariates are interpolated and affect PK", { + covs_locf <- list( + CRCL = new_covariate( + times = c(0, 48, 96), + value = c(3, 10, 3), + implementation = "locf" + ) + ) + covs_inter <- list( + CRCL = new_covariate( + times = c(0, 48, 96), + value = c(3, 10, 3), + implementation = "interpolate" + ) + ) + t_obs <- seq(0, 120, 2) + sim2_inter <- sim_ode( + mod_2cmt_timevar, + parameters = par_timevar, + regimen = reg_timevar, + covariates = covs_inter, + only_obs = TRUE, + t_obs = t_obs, + output_include = list(parameters = TRUE, covariates = TRUE, variables = TRUE) + ) + sim2_locf <- sim_ode( + mod_2cmt_timevar, + parameters = par_timevar, + regimen = reg_timevar, + covariates = covs_locf, + only_obs = TRUE, + t_obs = t_obs, + output_include = list(parameters = TRUE, covariates = TRUE, variables = TRUE) + ) + + ## Check covariate changing for inter, but not for locf + expect_equal( + round(sim2_inter$CRCL, 3)[1:10], + c(3, 3.292, 3.583, 3.875, 4.167, 4.458, 4.75, 5.042, 5.333, 5.625) + ) + expect_equal( + round(sim2_locf$CRCL, 3)[1:10], + rep(3, 10) + ) + + ## Check interpolated covariates actually affect PK parameters + expect_equal( + round(sim2_inter$CLi, 3)[1:10], + c(6, 6.292, 6.583, 6.875, 7.167, 7.458, 7.75, 8.042, 8.333, 8.625) + ) + expect_equal( + round(sim2_locf$CLi, 3)[1:10], + rep(6, 10) + ) + + ## Check interpolated covariates actually affect simulated conc + expect_equal( + round(sim2_inter$y, 3)[1:10], + c(0, 0.32, 0.611, 0.788, 1.205, 1.531, 1.708, 2.098, 2.379, 2.503) + ) + expect_equal( + round(sim2_locf$y, 3)[1:10], + c(0, 0.321, 0.617, 0.805, 1.239, 1.598, 1.813, 2.251, 2.598, 2.792) + ) + + ## Visual check: + # ggplot(sim2_inter, aes(x = t, y = y)) + + # geom_line() + + # geom_line(data = sim2_locf, colour = "blue") + + }) + +}) \ No newline at end of file diff --git a/tests/testthat/test_sim_lagtime.R b/tests/testthat/test_sim_lagtime.R deleted file mode 100644 index 77a8c5f7..00000000 --- a/tests/testthat/test_sim_lagtime.R +++ /dev/null @@ -1,26 +0,0 @@ -test_that("dose dump after lagtime in correct order in output data", { - skip_on_cran() - reg <- new_regimen(amt = 500, n = 4, interval = 12, type = 'oral') - pars <- list(CL = 5, V = 50, KA = 0.5, TLAG = 0.83) - dat <- sim_ode( - ode = mod_1cmt_oral_lagtime, - regimen = reg, - parameters = pars, - only_obs = FALSE - ) - ## Change after RXR-2394: time of TLAG not in dataset anymore unless requested by user in t_obs - ## Before change: expect_equal(round(dat[dat$t == 12.83 & dat$comp == 1,]$y, 1), c(1.2, 501.2)) - ## After change: - expect_equal(nrow(dat[dat$t == 12.83,]), 0) - ## When grid requested by user, lagtime should be visible - dat <- sim_ode( - ode = mod_1cmt_oral_lagtime, - regimen = reg, - parameters = pars, - t_obs = seq(0, 1, 0.01), - only_obs = TRUE - ) - tmp <- dat[dat$t >= 0.82 & dat$t <= 0.85, ] - expect_equal(tmp$t, c(0.82, 0.83, 0.84, 0.85)) - expect_equal(round(tmp$y, 2), c(0, 0, 0.05, 0.10)) -}) diff --git a/tests/testthat/test_t_init.R b/tests/testthat/test_t_init.R deleted file mode 100644 index f0e1a63a..00000000 --- a/tests/testthat/test_t_init.R +++ /dev/null @@ -1,37 +0,0 @@ -# Test t_init functionality - -## e.g. TDM before first dose: -## at t=-8, conc=10000 -## Use this as true value in the simulations -par <- list(CL = 7.67, V = 97.7, TDM_INIT = 500) -mod <- new_ode_model( - code = "CLi = CL; Vi = V; dAdt[1] = -(CLi/Vi)*A[1]; CONC = A[1]/Vi", - state_init = "A[1] = TDM_INIT * Vi", - parameters = par, - obs = list(cmt = 1, scale = "Vi"), - declare_variables = c("CONC", "CLi", "Vi"), - cpp_show_code = F -) -reg <- new_regimen(amt = 100000, times=c(0, 24), type="bolus") -s <- sim( - ode = mod, - parameters = par, - n_ind = 1, - regimen = reg, - only_obs = TRUE, - output_include = list(variables = TRUE), - t_init = 10 -) - -test_that("TDM before first dose is considered a true initial value", { - expect_equal(sum(is.na(s$y)), 0) - expect_equal(s$y[1], 500) - expect_equal(round(s$y[3]), 427) - expect_equal(round(s$y[13]),1157) -}) - -test_that("Variables are set (also in first row) when TDM before first dose", { - expect_equal(round(s$CONC[1:5], 1), c(500, 462.2, 427.3, 395.1, 365.3)) -}) - - diff --git a/tests/testthat/test_time_rounding_bug.R b/tests/testthat/test_time_rounding_bug.R deleted file mode 100644 index 018af52d..00000000 --- a/tests/testthat/test_time_rounding_bug.R +++ /dev/null @@ -1,15 +0,0 @@ -# rounding time should not produce NAs in sim -## time-rounding bug 20170804 - -test_that("No NAs related to rounding", { - # Uses models defined in setup.R - p <- list(CL = 5, V = 50, Q = 10, V2 = 150) - r1 <- new_regimen(amt = 100, times = c(0, 24, 36), type = "infusion") - dat1 <- sim_ode( - ode = mod_2cmt_iv, - parameters = p, - regimen = r1, - t_obs=seq(0, 150, length.out = 100) - ) - expect_equal(sum(is.na(dat1$y)), 0) -}) diff --git a/tests/testthat/test_timevar_cov.R b/tests/testthat/test_timevar_cov.R deleted file mode 100644 index 60866dbc..00000000 --- a/tests/testthat/test_timevar_cov.R +++ /dev/null @@ -1,118 +0,0 @@ -skip_on_cran() - -par <- list(CL = 3, V = 50, Q = 2.5, V2 = 70) -reg <- new_regimen( - amt = 250, - n = 60, - interval = 6, - type = 'infusion', - t_inf = 1 -) -mod <- new_ode_model( - code = " - dAdt[1] = -(Q/V)*A[1] + (Q/V2)*A[2] -(CLi/V)*A[1]; - dAdt[2] = -(Q/V2)*A[2] + (Q/V)*A[1]; - ", - pk_code = "CLi = CL + CRCL", - obs = list(cmt = 2, scale = "V"), - covariates = list(CRCL = new_covariate(5)), declare_variables = "CLi", - cpp = FALSE -) - -test_that("timevarying covariates handled", { - # CLi changes by several orders of magnitude after - # steady state is achieved, which should produce - # a new steady state that is much lower - covs <- list( - CRCL = new_covariate( - value = c(4.662744, 5.798767, 6.195943, 6.2, 600), - times = c(0, 18, 29, 209, 210), - implementation = "locf" - ) - ) - t_obs <- seq(0, 360, 0.1) - sim1 <- sim_ode( - mod, - parameters = par, - regimen = reg, - covariates = covs, - only_obs = TRUE, - t_obs = t_obs, - output_include = list(parameters = TRUE, covariates = TRUE) - ) - expect_equal(sim1$CRCL[sim1$t == 209], 6.2) - expect_equal(sim1$CRCL[sim1$t == 210], 600) - expect_true(sim1$y[sim1$t == 35 * 6] > 10 * sim1$y[sim1$t == 60 * 6]) -}) - -test_that("timevarying covariates are interpolated and affect PK", { - covs_locf <- list( - CRCL = new_covariate( - times = c(0, 48, 96), - value = c(3, 10, 3), - implementation = "locf" - ) - ) - covs_inter <- list( - CRCL = new_covariate( - times = c(0, 48, 96), - value = c(3, 10, 3), - implementation = "interpolate" - ) - ) - t_obs <- seq(0, 120, 2) - sim2_inter <- sim_ode( - mod, - parameters = par, - regimen = reg, - covariates = covs_inter, - only_obs = TRUE, - t_obs = t_obs, - output_include = list(parameters = TRUE, covariates = TRUE, variables = TRUE) - ) - sim2_locf <- sim_ode( - mod, - parameters = par, - regimen = reg, - covariates = covs_locf, - only_obs = TRUE, - t_obs = t_obs, - output_include = list(parameters = TRUE, covariates = TRUE, variables = TRUE) - ) - - ## Check covariate changing for inter, but not for locf - expect_equal( - round(sim2_inter$CRCL, 3)[1:10], - c(3, 3.292, 3.583, 3.875, 4.167, 4.458, 4.75, 5.042, 5.333, 5.625) - ) - expect_equal( - round(sim2_locf$CRCL, 3)[1:10], - rep(3, 10) - ) - - ## Check interpolated covariates actually affect PK parameters - expect_equal( - round(sim2_inter$CLi, 3)[1:10], - c(6, 6.292, 6.583, 6.875, 7.167, 7.458, 7.75, 8.042, 8.333, 8.625) - ) - expect_equal( - round(sim2_locf$CLi, 3)[1:10], - rep(6, 10) - ) - - ## Check interpolated covariates actually affect simulated conc - expect_equal( - round(sim2_inter$y, 3)[1:10], - c(0, 0.32, 0.611, 0.788, 1.205, 1.531, 1.708, 2.098, 2.379, 2.503) - ) - expect_equal( - round(sim2_locf$y, 3)[1:10], - c(0, 0.321, 0.617, 0.805, 1.239, 1.598, 1.813, 2.251, 2.598, 2.792) - ) - - ## Visual check: - # ggplot(sim2_inter, aes(x = t, y = y)) + - # geom_line() + - # geom_line(data = sim2_locf, colour = "blue") - -})