diff options
| -rw-r--r-- | R/mkinfit.R | 14 | ||||
| -rw-r--r-- | R/mkinmod.R | 2 | ||||
| -rw-r--r-- | R/mkinpredict.R | 13 | ||||
| -rw-r--r-- | R/plot.mixed.mmkin.R | 10 | ||||
| -rw-r--r-- | R/plot.mkinfit.R | 5 | ||||
| -rw-r--r-- | R/saem.R | 8 | ||||
| -rw-r--r-- | inst/rmarkdown/templates/hierarchical_kinetics/skeleton/skeleton.Rmd | 24 | ||||
| -rw-r--r-- | log/test.log | 29 | ||||
| -rw-r--r-- | man/mkinpredict.Rd | 4 | ||||
| -rw-r--r-- | tests/testthat/print_mmkin_sfo_sfo_dmta.txt | 0 | ||||
| -rw-r--r-- | tests/testthat/setup_script.R | 12 | ||||
| -rw-r--r-- | tests/testthat/test_compiled_symbols.R | 58 | ||||
| -rw-r--r-- | tests/testthat/test_dmta.R | 11 |
13 files changed, 145 insertions, 45 deletions
diff --git a/R/mkinfit.R b/R/mkinfit.R index b97bc7e2..c851fddb 100644 --- a/R/mkinfit.R +++ b/R/mkinfit.R @@ -501,10 +501,15 @@ mkinfit <- function(mkinmod, observed, } # Get native symbol before iterations info for speed + use_symbols = FALSE if (solution_type == "deSolve" & use_compiled[1] != FALSE) { - mkinmod[["symbols"]] <- deSolve::checkDLL(dllname = mkinmod$dll_info[["name"]], - func = "diffs", initfunc = "initpar", - jacfunc = NULL, nout = 0, outnames = NULL) + mkinmod[["symbols"]] <- try( + deSolve::checkDLL(dllname = mkinmod$dll_info[["name"]], + func = "diffs", initfunc = "initpar", + jacfunc = NULL, nout = 0, outnames = NULL)) + if (!inherits(mkinmod[["symbols"]], "try-error")) { + use_symbols = TRUE + } } # Get the error model and the algorithm for fitting @@ -616,8 +621,9 @@ mkinfit <- function(mkinmod, observed, odeini, outtimes, solution_type = solution_type, use_compiled = use_compiled, + use_symbols = use_symbols, method.ode = method.ode, - atol = atol, rtol = rtol, + atol = atol, rtol = rtol, ...) observed_index <- cbind(as.character(observed$time), as.character(observed$name)) diff --git a/R/mkinmod.R b/R/mkinmod.R index 29542b71..2f930adb 100644 --- a/R/mkinmod.R +++ b/R/mkinmod.R @@ -464,11 +464,11 @@ mkinmod <- function(..., use_of_ff = "max", name = NULL, silent = TRUE) if (!inherits(model$cf, "try-error")) { + if (!quiet) message("Temporary DLL for differentials generated and loaded") if (!is.null(dll_dir)) { model$dll_info <- inline::moveDLL(model$cf, name, dll_dir, unload = unload, overwrite = overwrite, verbose = !quiet) } - if (!quiet) message("Temporary DLL for differentials generated and loaded") model$dll_info <- inline::getDynLib(model$cf) } } diff --git a/R/mkinpredict.R b/R/mkinpredict.R index 957d5793..60456fb2 100644 --- a/R/mkinpredict.R +++ b/R/mkinpredict.R @@ -27,6 +27,8 @@ #' When using compiled code, only lsoda is supported. #' @param use_compiled If set to \code{FALSE}, no compiled version of the #' [mkinmod] model is used, even if is present. +#' @param use_symbols If set to \code{TRUE} (default), symbol info present in +#' the [mkinmod] object is used if available for accessing compiled code #' @param atol Absolute error tolerance, passed to the ode solver. #' @param rtol Absolute error tolerance, passed to the ode solver. #' @param maxsteps Maximum number of steps, passed to the ode solver. @@ -114,6 +116,7 @@ mkinpredict.mkinmod <- function(x, outtimes = seq(0, 120, by = 0.1), solution_type = "deSolve", use_compiled = "auto", + use_symbols = FALSE, method.ode = "lsoda", atol = 1e-8, rtol = 1e-10, maxsteps = 20000L, map_output = TRUE, na_stop = TRUE, @@ -169,12 +172,18 @@ mkinpredict.mkinmod <- function(x, } if (solution_type == "deSolve") { - if (!is.null(x$cf) & !is.null(x$symbols) & use_compiled[1] != FALSE) { + if (!is.null(x$cf) & use_compiled[1] != FALSE) { + + if (!is.null(x$symbols) & use_symbols) { + lsoda_func <- x$symbols + } else { + lsoda_func <- "diffs" + } out <- deSolve::lsoda( y = odeini, times = outtimes, - func = x$symbols, + func = lsoda_func, initfunc = "initpar", dllname = x$dll_info[["name"]], parms = odeparms[x$parms], # Order matters when using compiled models diff --git a/R/plot.mixed.mmkin.R b/R/plot.mixed.mmkin.R index 466e55fc..d6c3d0de 100644 --- a/R/plot.mixed.mmkin.R +++ b/R/plot.mixed.mmkin.R @@ -166,10 +166,12 @@ plot.mixed.mmkin <- function(x, } } - # Make sure degparms_pop is a matrix, columns corresponding to population curve(s) - if (is.null(dim(degparms_pop))) { - degparms_pop <- matrix(degparms_pop, ncol = 1, - dimnames = list(names(degparms_pop), "Population")) + if (pop_curves) { + # Make sure degparms_pop is a matrix, columns corresponding to population curve(s) + if (is.null(dim(degparms_pop))) { + degparms_pop <- matrix(degparms_pop, ncol = 1, + dimnames = list(names(degparms_pop), "Population")) + } } degparms_fixed <- fit_1$fixed$value diff --git a/R/plot.mkinfit.R b/R/plot.mkinfit.R index 9606507a..9f19213a 100644 --- a/R/plot.mkinfit.R +++ b/R/plot.mkinfit.R @@ -146,9 +146,8 @@ plot.mkinfit <- function(x, fit = x, func = "diffs", initfunc = "initpar", jacfunc = NULL, nout = 0, outnames = NULL) } - out <- try(mkinpredict(fit$mkinmod, odeparms, odeini, outtimes, - solution_type = solution_type, atol = fit$atol, rtol = fit$rtol), - silent = TRUE) + out <- mkinpredict(fit$mkinmod, odeparms, odeini, outtimes, + solution_type = solution_type, atol = fit$atol, rtol = fit$rtol) out <- as.data.frame(out) @@ -583,11 +583,15 @@ saemix_model <- function(object, solution_type = "auto", transform_fractions <- object[[1]]$transform_fractions # Get native symbol info for speed + use_symbols = FALSE if (solution_type == "deSolve" & !is.null(mkin_model$cf)) { - mkin_model$symbols <- deSolve::checkDLL( + mkin_model$symbols <- try(deSolve::checkDLL( dllname = mkin_model$dll_info[["name"]], func = "diffs", initfunc = "initpar", - jacfunc = NULL, nout = 0, outnames = NULL) + jacfunc = NULL, nout = 0, outnames = NULL)) + if (!inherits(mkin_model$symbols, "try-error")) { + use_symbols = TRUE + } } # Define the model function diff --git a/inst/rmarkdown/templates/hierarchical_kinetics/skeleton/skeleton.Rmd b/inst/rmarkdown/templates/hierarchical_kinetics/skeleton/skeleton.Rmd index 38a6bd20..cb328308 100644 --- a/inst/rmarkdown/templates/hierarchical_kinetics/skeleton/skeleton.Rmd +++ b/inst/rmarkdown/templates/hierarchical_kinetics/skeleton/skeleton.Rmd @@ -21,11 +21,18 @@ library(readxl) ```{r n_cores, cache = FALSE} n_cores <- detectCores() -if (Sys.info()["sysname"] == "Windows") { - cl <- makePSOCKcluster(n_cores) -} else { - cl <- makeForkCluster(n_cores) +# We need to start a new cluster after defining a compiled model that is +# saved as a DLL to the user directory, therefore we define a function +# This is used again after defining the pathway model +start_cluster <- function(n_cores) { + if (Sys.info()["sysname"] == "Windows") { + ret <- makePSOCKcluster(n_cores) + } else { + ret <- makeForkCluster(n_cores) + } + return(ret) } +cl <- start_cluster(n_cores) ``` \clearpage @@ -179,6 +186,10 @@ illparms(parent_best_pH_2) parms(parent_best_pH_2, ci = TRUE) |> kable(digits = 3) ``` +```{r} +stopCluster(cl) +``` + \clearpage # Pathway fits @@ -203,6 +214,7 @@ Separate evaluations of all datasets are performed with constant variance and using two-component error. ```{r path-1-sep, dependson = c("path-1-degmod", "ds")} +cl <- start_cluster(n_cores) sforb_sep_const <- mmkin(list(sforb_path = m_sforb_sfo2), ds, cluster = cl, quiet = TRUE) sforb_sep_tc <- update(sforb_sep_const, error_model = "tc") @@ -268,6 +280,10 @@ plot(path_1_refined) parms(path_1_refined, ci = TRUE) |> kable(digits = 3) ``` +```{r} +stopCluster(cl) +``` + \clearpage # Appendix diff --git a/log/test.log b/log/test.log index c89cf02e..fa808f22 100644 --- a/log/test.log +++ b/log/test.log @@ -4,20 +4,21 @@ ✔ | 5 | Analytical solutions for coupled models [1.6s] ✔ | 5 | Calculation of Akaike weights ✔ | 3 | Export dataset for reading into CAKE +✔ | 7 | Use of precompiled symbols in mkinpredict [3.3s] ✔ | 12 | Confidence intervals and p-values [0.4s] -✔ | 1 12 | Dimethenamid data from 2018 [12.2s] +✔ | 1 12 | Dimethenamid data from 2018 [13.4s] ──────────────────────────────────────────────────────────────────────────────── -Skip ('test_dmta.R:99'): Different backends get consistent results for SFO-SFO3+, dimethenamid data +Skip ('test_dmta.R:88'): 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 [2.4s] +✔ | 14 | Error model fitting [2.5s] ✔ | 5 | Time step normalisation ✔ | 4 | Calculation of FOCUS chi2 error levels [0.3s] -✔ | 14 | Results for FOCUS D established in expertise for UBA (Ranke 2014) [0.4s] -✔ | 4 | Test fitting the decline of metabolites from their maximum [0.2s] +✔ | 14 | Results for FOCUS D established in expertise for UBA (Ranke 2014) [0.5s] +✔ | 4 | Test fitting the decline of metabolites from their maximum [0.3s] ✔ | 1 | Fitting the logistic model [0.1s] -✔ | 10 | Batch fitting and diagnosing hierarchical kinetic models [19.5s] -✔ | 1 11 | Nonlinear mixed-effects models [5.8s] +✔ | 10 | Batch fitting and diagnosing hierarchical kinetic models [19.8s] +✔ | 1 11 | Nonlinear mixed-effects models [6.1s] ──────────────────────────────────────────────────────────────────────────────── Skip ('test_mixed.R:78'): saemix results are reproducible for biphasic fits Reason: Fitting with saemix takes around 10 minutes when using deSolve @@ -27,12 +28,12 @@ Reason: Fitting with saemix takes around 10 minutes when using deSolve ✔ | 3 | mkinfit features [0.5s] ✔ | 8 | mkinmod model generation and printing ✔ | 3 | Model predictions with mkinpredict [0.1s] -✔ | 12 | Multistart method for saem.mmkin models [22.0s] -✔ | 16 | Evaluations according to 2015 NAFTA guidance [1.6s] -✔ | 9 | Nonlinear mixed-effects models with nlme [3.9s] -✔ | 15 | Plotting [5.1s] +✔ | 12 | Multistart method for saem.mmkin models [21.7s] +✔ | 16 | Evaluations according to 2015 NAFTA guidance [1.5s] +✔ | 9 | Nonlinear mixed-effects models with nlme [3.8s] +✔ | 15 | Plotting [4.7s] ✔ | 4 | Residuals extracted from mkinfit models -✔ | 1 36 | saemix parent models [31.4s] +✔ | 1 36 | saemix parent models [32.0s] ──────────────────────────────────────────────────────────────────────────────── Skip ('test_saemix_parent.R:143'): We can also use mkin solution methods for saem Reason: This still takes almost 2.5 minutes although we do not solve ODEs @@ -47,11 +48,11 @@ Reason: This still takes almost 2.5 minutes although we do not solve ODEs ✔ | 4 | Calculation of maximum time weighted average concentrations (TWAs) [0.7s] ══ Results ═════════════════════════════════════════════════════════════════════ -Duration: 115.0 s +Duration: 120.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) • This still takes almost 2.5 minutes although we do not solve ODEs (1) -[ FAIL 0 | WARN 0 | SKIP 3 | PASS 270 ] +[ FAIL 0 | WARN 0 | SKIP 3 | PASS 277 ] diff --git a/man/mkinpredict.Rd b/man/mkinpredict.Rd index ad749fa4..792d0e47 100644 --- a/man/mkinpredict.Rd +++ b/man/mkinpredict.Rd @@ -15,6 +15,7 @@ mkinpredict(x, odeparms, odeini, outtimes, ...) outtimes = seq(0, 120, by = 0.1), solution_type = "deSolve", use_compiled = "auto", + use_symbols = FALSE, method.ode = "lsoda", atol = 1e-08, rtol = 1e-10, @@ -67,6 +68,9 @@ parent compound.} \item{use_compiled}{If set to \code{FALSE}, no compiled version of the \link{mkinmod} model is used, even if is present.} +\item{use_symbols}{If set to \code{TRUE} (default), symbol info present in +the \link{mkinmod} object is used if available for accessing compiled code} + \item{method.ode}{The solution method passed via \link{mkinpredict} to \link{ode}] in case the solution type is "deSolve" and we are not using compiled code. When using compiled code, only lsoda is supported.} diff --git a/tests/testthat/print_mmkin_sfo_sfo_dmta.txt b/tests/testthat/print_mmkin_sfo_sfo_dmta.txt new file mode 100644 index 00000000..e69de29b --- /dev/null +++ b/tests/testthat/print_mmkin_sfo_sfo_dmta.txt diff --git a/tests/testthat/setup_script.R b/tests/testthat/setup_script.R index 4e2f76ab..b2147fbe 100644 --- a/tests/testthat/setup_script.R +++ b/tests/testthat/setup_script.R @@ -110,3 +110,15 @@ dfop_saem_1 <- saem(mmkin_dfop_1, quiet = TRUE, transformations = "mkin", no_random_effect = c("parent_0", "g_qlogis")) parallel::stopCluster(cl) + +# Preprocess dimethenamid data +dmta_ds <- lapply(1:7, function(i) { + ds_i <- dimethenamid_2018$ds[[i]]$data + ds_i[ds_i$name == "DMTAP", "name"] <- "DMTA" + ds_i$time <- ds_i$time * dimethenamid_2018$f_time_norm[i] + ds_i +}) +names(dmta_ds) <- sapply(dimethenamid_2018$ds, function(ds) ds$title) +dmta_ds[["Elliot"]] <- rbind(dmta_ds[["Elliot 1"]], dmta_ds[["Elliot 2"]]) +dmta_ds[["Elliot 1"]] <- dmta_ds[["Elliot 2"]] <- NULL + diff --git a/tests/testthat/test_compiled_symbols.R b/tests/testthat/test_compiled_symbols.R new file mode 100644 index 00000000..2828951d --- /dev/null +++ b/tests/testthat/test_compiled_symbols.R @@ -0,0 +1,58 @@ +context("Use of precompiled symbols in mkinpredict") + +test_that("We can safely use compiled code", { + + # Generate temporary DLL + sfo_sfo_tmp <- mkinmod(DMTA = mkinsub("SFO", to = "M23"), + M23 = mkinsub("SFO")) + + # Generate temporary DLL and move to user specified location + if (!dir.exists("test_dlls")) dir.create("test_dlls") + sfo_sfo_dll <- mkinmod(DMTA = mkinsub("SFO", to = "M23"), + M23 = mkinsub("SFO"), + dll_dir = "test_dlls", + name = "sfo_sfo", + unload = TRUE, overwrite = TRUE + ) + + if (Sys.info()["sysname"] != "Windows") { + # mclapply using forks + expect_known_output( + mmkin(list(sfo_sfo_dll), dmta_ds, cores = n_cores, quiet = TRUE), + "print_mmkin_sfo_sfo_dmta.txt" + ) + + # cluster describing itself as socket cluster + cl_fork <- parallel::makeForkCluster(n_cores) + expect_known_output( + mmkin(list(sfo_sfo_tmp), dmta_ds, cluster = cl_fork, quiet = TRUE), + "print_mmkin_sfo_sfo_dmta.txt" + ) + expect_known_output( + mmkin(list(sfo_sfo_dll), dmta_ds, cluster = cl_fork, quiet = TRUE), + "print_mmkin_sfo_sfo_dmta.txt" + ) + parallel::stopCluster(cl_fork) + } + + # PSOCK cluster + cl_psock <- parallel::makePSOCKcluster(n_cores) + expect_known_output( + mmkin(list(sfo_sfo_tmp), dmta_ds, cluster = cl_psock, quiet = TRUE), + "print_mmkin_sfo_sfo_dmta.txt" + ) + expect_known_output( + mmkin(list(sfo_sfo_dll), dmta_ds, cluster = cl_psock, quiet = TRUE), + "print_mmkin_sfo_sfo_dmta.txt" + ) + parallel::stopCluster(cl_psock) + + # Clean up + if (Sys.info()["sysname"] != "Windows") { + expect_true(file.remove("test_dlls/sfo_sfo.dll")) + } else { + expect_true(file.remove("test_dlls/sfo_sfo.so")) + } + expect_true(file.remove("test_dlls")) +}) + diff --git a/tests/testthat/test_dmta.R b/tests/testthat/test_dmta.R index 825c6e80..cd6de341 100644 --- a/tests/testthat/test_dmta.R +++ b/tests/testthat/test_dmta.R @@ -1,16 +1,5 @@ context("Dimethenamid data from 2018") -# Data -dmta_ds <- lapply(1:7, function(i) { - ds_i <- dimethenamid_2018$ds[[i]]$data - ds_i[ds_i$name == "DMTAP", "name"] <- "DMTA" - ds_i$time <- ds_i$time * dimethenamid_2018$f_time_norm[i] - ds_i -}) -names(dmta_ds) <- sapply(dimethenamid_2018$ds, function(ds) ds$title) -dmta_ds[["Elliot"]] <- rbind(dmta_ds[["Elliot 1"]], dmta_ds[["Elliot 2"]]) -dmta_ds[["Elliot 1"]] <- dmta_ds[["Elliot 2"]] <- NULL - test_that("Different backends get consistent results for DFOP tc, dimethenamid data", { skip_on_cran() # Time constraints |