From 2b331b5f9420943ebe26b60361578038ec560f88 Mon Sep 17 00:00:00 2001 From: Johannes Ranke Date: Tue, 1 Mar 2022 11:23:49 +0100 Subject: Add some more tests --- R/nlme.mmkin.R | 3 +- man/saem.Rd | 1 + print_nlmixr_saem_biphasic.txt | 35 -------------- print_nlmixr_saem_biphasic_tc.txt | 39 ---------------- test.log | 28 +++++++----- tests/testthat/dimethenamid_2018.txt | 24 ++++++++++ tests/testthat/print_nlmixr_saem_biphasic.txt | 35 ++++++++++++++ tests/testthat/print_nlmixr_saem_biphasic_tc.txt | 40 ++++++++++++++++ tests/testthat/setup_script.R | 8 ++-- tests/testthat/summary_nlmixr_saem_biphasic.txt | 58 ++++++++++++------------ tests/testthat/test_dmta.R | 5 +- tests/testthat/test_mkinds.R | 2 + 12 files changed, 157 insertions(+), 121 deletions(-) delete mode 100644 print_nlmixr_saem_biphasic.txt delete mode 100644 print_nlmixr_saem_biphasic_tc.txt create mode 100644 tests/testthat/dimethenamid_2018.txt create mode 100644 tests/testthat/print_nlmixr_saem_biphasic.txt create mode 100644 tests/testthat/print_nlmixr_saem_biphasic_tc.txt diff --git a/R/nlme.mmkin.R b/R/nlme.mmkin.R index 79139d12..b6b6c830 100644 --- a/R/nlme.mmkin.R +++ b/R/nlme.mmkin.R @@ -193,7 +193,8 @@ nlme.mmkin <- function(model, data = "auto", val$mkinmod <- model[[1]]$mkinmod val$data <- thisCall[["data"]] val$mmkin <- model - val$mean_dp_start <- ifelse(is.list(start), start$fixed, start) + if (is.list(start)) val$mean_dp_start <- start$fixed + else val$mean_dp_start <- start val$transform_rates <- model[[1]]$transform_rates val$transform_fractions <- model[[1]]$transform_fractions val$solution_type <- model[[1]]$solution_type diff --git a/man/saem.Rd b/man/saem.Rd index 00e9aeda..0c066dd2 100644 --- a/man/saem.Rd +++ b/man/saem.Rd @@ -34,6 +34,7 @@ saemix_model( transformations = c("mkin", "saemix"), degparms_start = numeric(), test_log_parms = FALSE, + conf.level = 0.6, verbose = FALSE, ... ) diff --git a/print_nlmixr_saem_biphasic.txt b/print_nlmixr_saem_biphasic.txt deleted file mode 100644 index 5cf80f0f..00000000 --- a/print_nlmixr_saem_biphasic.txt +++ /dev/null @@ -1,35 +0,0 @@ -Kinetic nonlinear mixed-effects model fit by saem using nlmixr -Structural model: -d_parent/dt = - ((k1 * g * exp(-k1 * time) + k2 * (1 - g) * exp(-k2 * - time)) / (g * exp(-k1 * time) + (1 - g) * exp(-k2 * time))) - * parent -d_m1/dt = + f_parent_to_m1 * ((k1 * g * exp(-k1 * time) + k2 * (1 - g) - * exp(-k2 * time)) / (g * exp(-k1 * time) + (1 - g) * - exp(-k2 * time))) * parent - k_m1 * m1 - -Data: -507 observations of 2 variable(s) grouped in 15 datasets - -Likelihood: - AIC BIC logLik - 2510 2569 -1241 - -Fitted parameters: - Est. SE %RSE Back-transformed(95%CI) BSV(CV% or SD) -parent_0 100 0.444 0.442 100 (99.6, 101) 0.461 -log_k_m1 -5.36 0.0804 1.5 0.00472 (0.00403, 0.00552) 21.5 -f_parent_qlogis -0.0909 0.0868 95.5 0.477 (0.435, 0.52) 0.308 -log_k1 -2.78 0.0977 3.51 0.0618 (0.0511, 0.0749) 30.1 -log_k2 -4.51 0.08 1.77 0.011 (0.00941, 0.0129) 16.6 -g_qlogis -0.121 0.17 141 0.47 (0.388, 0.553) 0.373 -sigma_parent 3.72 3.72 -sigma_m1 1.44 1.44 - Shrink(SD)% -parent_0 69.2%> -log_k_m1 19.6%< -f_parent_qlogis -0.136%> -log_k1 6.71%< -log_k2 23.7%= -g_qlogis 19.5%< -sigma_parent -sigma_m1 diff --git a/print_nlmixr_saem_biphasic_tc.txt b/print_nlmixr_saem_biphasic_tc.txt deleted file mode 100644 index ad59005d..00000000 --- a/print_nlmixr_saem_biphasic_tc.txt +++ /dev/null @@ -1,39 +0,0 @@ -Kinetic nonlinear mixed-effects model fit by saem using nlmixr -Structural model: -d_parent/dt = - ((k1 * g * exp(-k1 * time) + k2 * (1 - g) * exp(-k2 * - time)) / (g * exp(-k1 * time) + (1 - g) * exp(-k2 * time))) - * parent -d_m1/dt = + f_parent_to_m1 * ((k1 * g * exp(-k1 * time) + k2 * (1 - g) - * exp(-k2 * time)) / (g * exp(-k1 * time) + (1 - g) * - exp(-k2 * time))) * parent - k_m1 * m1 - -Data: -507 observations of 2 variable(s) grouped in 15 datasets - -Likelihood: - AIC BIC logLik - 2400 2468 -1184 - -Fitted parameters: - Est. SE %RSE Back-transformed(95%CI) BSV(CV% or SD) -parent_0 100 0.685 0.682 100 (99, 102) 0.278 -log_k_m1 -5.35 0.0812 1.52 0.00477 (0.00407, 0.00559) 13.7 -f_parent_qlogis -0.091 0.096 106 0.477 (0.431, 0.524) 0.320 -log_k1 -2.76 0.0969 3.5 0.063 (0.0521, 0.0762) 29.9 -log_k2 -4.49 0.0658 1.47 0.0113 (0.00991, 0.0128) 18.6 -g_qlogis -0.167 0.142 85.3 0.458 (0.39, 0.528) 0.339 -sigma_low_parent 1.04 1.04 -rsd_high_parent 0.0531 0.0531 -sigma_low_m1 0.819 0.819 -rsd_high_m1 0.058 0.058 - Shrink(SD)% -parent_0 87.7%> -log_k_m1 44.8%> -f_parent_qlogis -0.375%> -log_k1 6.76%< -log_k2 10.4%< -g_qlogis 20.3%= -sigma_low_parent -rsd_high_parent -sigma_low_m1 -rsd_high_m1 diff --git a/test.log b/test.log index 5c8f1802..b044e6a0 100644 --- a/test.log +++ b/test.log @@ -3,45 +3,51 @@ Loading required package: parallel ℹ Testing mkin ✔ | F W S OK | Context ✔ | 5 | AIC calculation -✔ | 5 | Analytical solutions for coupled models [3.4s] +✔ | 5 | Analytical solutions for coupled models [4.3s] ✔ | 5 | Calculation of Akaike weights ✔ | 2 | Export dataset for reading into CAKE -✔ | 12 | Confidence intervals and p-values [1.0s] +✔ | 12 | Confidence intervals and p-values [1.1s] ⠋ | 1 | Dimethenamid data from 2018 -✔ | 24 | Dimethenamid data from 2018 [37.1s] -✔ | 14 | Error model fitting [6.8s] +✔ | 1 25 | Dimethenamid data from 2018 [49.7s] +──────────────────────────────────────────────────────────────────────────────── +Skip (test_dmta.R:147:3): Different backends get consistent results for SFO-SFO3+, dimethenamid data +Reason: Fitting this ODE model with saemix takes about 15 minutes on my system +──────────────────────────────────────────────────────────────────────────────── +✔ | 14 | Error model fitting [7.0s] ✔ | 5 | Time step normalisation ✔ | 4 | Calculation of FOCUS chi2 error levels [0.8s] ✔ | 14 | Results for FOCUS D established in expertise for UBA (Ranke 2014) [1.1s] ✔ | 4 | Test fitting the decline of metabolites from their maximum [0.5s] ✔ | 1 | Fitting the logistic model [0.3s] -✔ | 1 12 | Nonlinear mixed-effects models [0.2s] +⠴ | 6 | Nonlinear mixed-effects models +✔ | 1 15 | Nonlinear mixed-effects models [2.4s] ──────────────────────────────────────────────────────────────────────────────── Skip (test_mixed.R:67:3): saemix results are reproducible for biphasic fits Reason: Fitting with saemix takes around 10 minutes when using deSolve ──────────────────────────────────────────────────────────────────────────────── -✔ | 2 | Test dataset classes mkinds and mkindsg +✔ | 3 | Test dataset classes mkinds and mkindsg ✔ | 10 | Special cases of mkinfit calls [0.6s] ✔ | 1 | mkinfit features [0.5s] ✔ | 8 | mkinmod model generation and printing [0.2s] ✔ | 3 | Model predictions with mkinpredict [0.4s] ✔ | 16 | Evaluations according to 2015 NAFTA guidance [2.1s] -✔ | 9 | Nonlinear mixed-effects models with nlme [9.0s] +✔ | 9 | Nonlinear mixed-effects models with nlme [9.2s] ✔ | 16 | Plotting [1.5s] ✔ | 4 | Residuals extracted from mkinfit models -✔ | 23 | saemix_parent [29.0s] +✔ | 23 | saemix_parent [30.4s] ✔ | 2 | Complex test case from Schaefer et al. (2007) Piacenza paper [1.6s] ✔ | 7 | Fitting the SFORB model [4.5s] ✔ | 1 | Summaries of old mkinfit objects ✔ | 4 | Summary [0.1s] -✔ | 4 | Results for synthetic data established in expertise for UBA (Ranke 2014) [2.5s] +✔ | 4 | Results for synthetic data established in expertise for UBA (Ranke 2014) [2.6s] ✔ | 9 | Hypothesis tests [9.4s] ✔ | 4 | Calculation of maximum time weighted average concentrations (TWAs) [2.2s] ══ Results ═════════════════════════════════════════════════════════════════════ -Duration: 115.2 s +Duration: 133.2 s ── Skipped tests ────────────────────────────────────────────────────────────── +• Fitting this ODE model with saemix takes about 15 minutes on my system (1) • Fitting with saemix takes around 10 minutes when using deSolve (1) -[ FAIL 0 | WARN 0 | SKIP 1 | PASS 230 ] +[ FAIL 0 | WARN 0 | SKIP 2 | PASS 235 ] diff --git a/tests/testthat/dimethenamid_2018.txt b/tests/testthat/dimethenamid_2018.txt new file mode 100644 index 00000000..513acff9 --- /dev/null +++ b/tests/testthat/dimethenamid_2018.txt @@ -0,0 +1,24 @@ + holding 7 mkinds objects +Title $title: Aerobic soil degradation data on dimethenamid-P from the EU assessment in 2018 +Occurrence of observed compounds $observed_n: +DMTAP M23 M27 M31 DMTA + 3 7 7 7 4 +Time normalisation factors $f_time_norm: +[1] 1.0000000 0.9706477 1.2284784 1.2284784 0.6233856 0.7678922 0.6733938 +Meta information $meta: + study usda_soil_type study_moisture_ref_type rel_moisture +Calke Unsworth 2014 Sandy loam pF2 1.00 +Borstel Staudenmaier 2009 Sand pF1 0.50 +Elliot 1 Wendt 1997 Clay loam pF2.5 0.75 +Elliot 2 Wendt 1997 Clay loam pF2.5 0.75 +Flaach König 1996 Sandy clay loam pF1 0.40 +BBA 2.2 König 1995 Loamy sand pF1 0.40 +BBA 2.3 König 1995 Sandy loam pF1 0.40 + study_ref_moisture temperature +Calke NA 20 +Borstel 23.00 20 +Elliot 1 33.37 23 +Elliot 2 33.37 23 +Flaach NA 20 +BBA 2.2 NA 20 +BBA 2.3 NA 20 diff --git a/tests/testthat/print_nlmixr_saem_biphasic.txt b/tests/testthat/print_nlmixr_saem_biphasic.txt new file mode 100644 index 00000000..5cf80f0f --- /dev/null +++ b/tests/testthat/print_nlmixr_saem_biphasic.txt @@ -0,0 +1,35 @@ +Kinetic nonlinear mixed-effects model fit by saem using nlmixr +Structural model: +d_parent/dt = - ((k1 * g * exp(-k1 * time) + k2 * (1 - g) * exp(-k2 * + time)) / (g * exp(-k1 * time) + (1 - g) * exp(-k2 * time))) + * parent +d_m1/dt = + f_parent_to_m1 * ((k1 * g * exp(-k1 * time) + k2 * (1 - g) + * exp(-k2 * time)) / (g * exp(-k1 * time) + (1 - g) * + exp(-k2 * time))) * parent - k_m1 * m1 + +Data: +507 observations of 2 variable(s) grouped in 15 datasets + +Likelihood: + AIC BIC logLik + 2510 2569 -1241 + +Fitted parameters: + Est. SE %RSE Back-transformed(95%CI) BSV(CV% or SD) +parent_0 100 0.444 0.442 100 (99.6, 101) 0.461 +log_k_m1 -5.36 0.0804 1.5 0.00472 (0.00403, 0.00552) 21.5 +f_parent_qlogis -0.0909 0.0868 95.5 0.477 (0.435, 0.52) 0.308 +log_k1 -2.78 0.0977 3.51 0.0618 (0.0511, 0.0749) 30.1 +log_k2 -4.51 0.08 1.77 0.011 (0.00941, 0.0129) 16.6 +g_qlogis -0.121 0.17 141 0.47 (0.388, 0.553) 0.373 +sigma_parent 3.72 3.72 +sigma_m1 1.44 1.44 + Shrink(SD)% +parent_0 69.2%> +log_k_m1 19.6%< +f_parent_qlogis -0.136%> +log_k1 6.71%< +log_k2 23.7%= +g_qlogis 19.5%< +sigma_parent +sigma_m1 diff --git a/tests/testthat/print_nlmixr_saem_biphasic_tc.txt b/tests/testthat/print_nlmixr_saem_biphasic_tc.txt new file mode 100644 index 00000000..2f940bfe --- /dev/null +++ b/tests/testthat/print_nlmixr_saem_biphasic_tc.txt @@ -0,0 +1,40 @@ +Kinetic nonlinear mixed-effects model fit by saem using nlmixr +Structural model: +d_parent/dt = - ((k1 * g * exp(-k1 * time) + k2 * (1 - g) * exp(-k2 * + time)) / (g * exp(-k1 * time) + (1 - g) * exp(-k2 * time))) + * parent +d_m1/dt = + f_parent_to_m1 * ((k1 * g * exp(-k1 * time) + k2 * (1 - g) + * exp(-k2 * time)) / (g * exp(-k1 * time) + (1 - g) * + exp(-k2 * time))) * parent - k_m1 * m1 + +Data: +507 observations of 2 variable(s) grouped in 15 datasets + +Likelihood: + + AIC BIC logLik + 2400 2468 -1184 + +Fitted parameters: + Est. SE %RSE Back-transformed(95%CI) BSV(CV% or SD) +parent_0 100 0.685 0.682 100 (99, 102) 0.278 +log_k_m1 -5.35 0.0812 1.52 0.00477 (0.00407, 0.00559) 13.7 +f_parent_qlogis -0.091 0.096 106 0.477 (0.431, 0.524) 0.320 +log_k1 -2.76 0.0969 3.5 0.063 (0.0521, 0.0762) 29.9 +log_k2 -4.49 0.0658 1.47 0.0113 (0.00991, 0.0128) 18.6 +g_qlogis -0.167 0.142 85.3 0.458 (0.39, 0.528) 0.339 +sigma_low_parent 1.04 1.04 +rsd_high_parent 0.0531 0.0531 +sigma_low_m1 0.819 0.819 +rsd_high_m1 0.058 0.058 + Shrink(SD)% +parent_0 87.7%> +log_k_m1 44.8%> +f_parent_qlogis -0.375%> +log_k1 6.76%< +log_k2 10.4%< +g_qlogis 20.3%= +sigma_low_parent +rsd_high_parent +sigma_low_m1 +rsd_high_m1 diff --git a/tests/testthat/setup_script.R b/tests/testthat/setup_script.R index bfa70005..4b99ad67 100644 --- a/tests/testthat/setup_script.R +++ b/tests/testthat/setup_script.R @@ -204,10 +204,10 @@ saem_biphasic_m <- saem(mmkin_biphasic, transformations = "mkin", quiet = TRUE) saem_biphasic_s <- saem(mmkin_biphasic, transformations = "saemix", quiet = TRUE) # nlmixr saem -tmp <- capture_output(nlmixr_saem_biphasic <- nlmixr(mmkin_biphasic, est = "saem", - control = nlmixr::saemControl(nBurn = 300, nEm = 100, nmc = 9, print = 0))) -tmp <- capture_output(nlmixr_saem_biphasic_tc <- nlmixr(mmkin_biphasic_tc, est = "saem", - control = nlmixr::saemControl(nBurn = 300, nEm = 100, nmc = 9, print = 0))) +tmp <- suppressMessages(capture.output(nlmixr_saem_biphasic <- nlmixr(mmkin_biphasic, est = "saem", + control = nlmixr::saemControl(nBurn = 300, nEm = 100, nmc = 9, print = 0)))) +tmp <- suppressMessages(capture.output(nlmixr_saem_biphasic_tc <- nlmixr(mmkin_biphasic_tc, est = "saem", + control = nlmixr::saemControl(nBurn = 300, nEm = 100, nmc = 9, print = 0)))) # The FOCEI fit takes too long... #tmp <- capture_output(nlmixr_focei_biphasic <- nlmixr(mmkin_biphasic, est = "focei", # control = nlmixr::foceiControl(print = 0))) diff --git a/tests/testthat/summary_nlmixr_saem_biphasic.txt b/tests/testthat/summary_nlmixr_saem_biphasic.txt index c5c1fe78..f1cdc47b 100644 --- a/tests/testthat/summary_nlmixr_saem_biphasic.txt +++ b/tests/testthat/summary_nlmixr_saem_biphasic.txt @@ -13,7 +13,7 @@ d_m1/dt = + f_parent_to_m1 * ((k1 * g * exp(-k1 * time) + k2 * (1 - g) exp(-k2 * time))) * parent - k_m1 * m1 Data: -508 observations of 2 variable(s) grouped in 15 datasets +507 observations of 2 variable(s) grouped in 15 datasets Degradation model predictions using RxODE @@ -23,9 +23,9 @@ Variance model: Constant variance Mean of starting values for individual parameters: parent_0 log_k_m1 f_parent_qlogis log_k1 log_k2 - 100.74 -5.43 -0.05 -2.85 -4.57 + 100.67 -5.38 -0.09 -2.74 -4.51 g_qlogis - -0.15 + -0.18 Mean of starting values for error model parameters: sigma_parent sigma_m1 @@ -38,39 +38,39 @@ Results: Likelihood calculated by gauss3_1.6 AIC BIC logLik - 2465 2525 -1219 + 2510 2569 -1241 Optimised parameters: est. lower upper -parent_0 100.54 99.7 101.40 -log_k_m1 -5.44 -5.6 -5.24 -f_parent_qlogis -0.06 -0.2 0.05 -log_k1 -2.81 -3.0 -2.65 -log_k2 -4.62 -4.8 -4.41 -g_qlogis -0.14 -0.4 0.17 +parent_0 100.49 99.6 101.37 +log_k_m1 -5.36 -5.5 -5.20 +f_parent_qlogis -0.09 -0.3 0.08 +log_k1 -2.78 -3.0 -2.59 +log_k2 -4.51 -4.7 -4.35 +g_qlogis -0.12 -0.5 0.21 Correlation: - prnt_0 lg_k_1 f_prn_ log_k1 log_k2 -log_k_m1 -0.224 -f_parent_qlogis -0.401 0.335 -log_k1 0.119 0.083 0.019 -log_k2 -0.019 0.237 0.142 0.395 -g_qlogis 0.083 -0.282 -0.218 -0.532 -0.544 + pr_0 l__1 f_p_ lg_1 lg_2 +log_k_m1 -0.3 +f_parent_qlogis -0.2 0.3 +log_k1 0.1 0.1 0.0 +log_k2 0.0 0.3 0.1 0.4 +g_qlogis 0.0 -0.3 -0.1 -0.5 -0.7 Random effects (omega): eta.parent_0 eta.log_k_m1 eta.f_parent_qlogis eta.log_k1 -eta.parent_0 0.05 0.0 0.00 0.00 -eta.log_k_m1 0.00 0.1 0.00 0.00 -eta.f_parent_qlogis 0.00 0.0 0.03 0.00 -eta.log_k1 0.00 0.0 0.00 0.05 -eta.log_k2 0.00 0.0 0.00 0.00 -eta.g_qlogis 0.00 0.0 0.00 0.00 +eta.parent_0 0.2 0.00 0.00 0.00 +eta.log_k_m1 0.0 0.05 0.00 0.00 +eta.f_parent_qlogis 0.0 0.00 0.09 0.00 +eta.log_k1 0.0 0.00 0.00 0.09 +eta.log_k2 0.0 0.00 0.00 0.00 +eta.g_qlogis 0.0 0.00 0.00 0.00 eta.log_k2 eta.g_qlogis eta.parent_0 0.00 0.0 eta.log_k_m1 0.00 0.0 eta.f_parent_qlogis 0.00 0.0 eta.log_k1 0.00 0.0 -eta.log_k2 0.08 0.0 +eta.log_k2 0.03 0.0 eta.g_qlogis 0.00 0.1 Variance model: @@ -80,11 +80,11 @@ sigma_parent sigma_m1 Backtransformed parameters: est. lower upper parent_0 1e+02 1e+02 1e+02 -k_m1 4e-03 4e-03 5e-03 -f_parent_to_m1 5e-01 5e-01 5e-01 +k_m1 5e-03 4e-03 6e-03 +f_parent_to_m1 5e-01 4e-01 5e-01 k1 6e-02 5e-02 7e-02 -k2 1e-02 8e-03 1e-02 -g 5e-01 4e-01 5e-01 +k2 1e-02 9e-03 1e-02 +g 5e-01 4e-01 6e-01 Resulting formation fractions: ff @@ -93,5 +93,5 @@ parent_sink 0.5 Estimated disappearance times: DT50 DT90 DT50back DT50_k1 DT50_k2 -parent 27 170 51 12 70 -m1 160 532 NA NA NA +parent 25 152 46 11 63 +m1 147 488 NA NA NA diff --git a/tests/testthat/test_dmta.R b/tests/testthat/test_dmta.R index 22fa9d95..7b130999 100644 --- a/tests/testthat/test_dmta.R +++ b/tests/testthat/test_dmta.R @@ -139,8 +139,9 @@ dmta_sfo_sfo3p_tc <- mmkin(list("SFO-SFO3+" = sfo_sfo3p), test_that("Different backends get consistent results for SFO-SFO3+, dimethenamid data", { - nlme_sfo_sfo3p_tc <- nlme(dmta_sfo_sfo3p_tc, - start = mean_degparms(dmta_sfo_sfo3p_tc, test_log_parms = TRUE)) + expect_warning(nlme_sfo_sfo3p_tc <- nlme(dmta_sfo_sfo3p_tc, + start = mean_degparms(dmta_sfo_sfo3p_tc, test_log_parms = TRUE)), + "Iteration 5, LME step.*not converge") ints_nlme_mets <- intervals(nlme_sfo_sfo3p_tc, which = "fixed") skip("Fitting this ODE model with saemix takes about 15 minutes on my system") diff --git a/tests/testthat/test_mkinds.R b/tests/testthat/test_mkinds.R index 5ce3619b..c10a3f5b 100644 --- a/tests/testthat/test_mkinds.R +++ b/tests/testthat/test_mkinds.R @@ -8,4 +8,6 @@ test_that("An mkinds object can be created and printed", { test_that("An mkindsg object can be created and printed", { testdata_group <- mkindsg$new("Experimental X", experimental_data_for_UBA_2019[6:10]) expect_known_output(print(testdata_group), "experimental_data_for_UBA_2019_mkindsg.txt") + + expect_known_output(print(dimethenamid_2018), "dimethenamid_2018.txt") }) -- cgit v1.2.1