From 1195dfc8bdbf7c131d6c6ec30fedbbe746af1bee Mon Sep 17 00:00:00 2001 From: Johannes Ranke Date: Wed, 6 May 2020 21:33:12 +0200 Subject: Change implementation of analytical solutions Preparing for symbolic solutions for more than one compound --- DESCRIPTION | 2 +- R/mkinfit.R | 6 +- R/mkinmod.R | 40 +++++++++++ R/mkinpredict.R | 152 ++++++++++++++++------------------------ build.log | 1 - check.log | 3 +- man/mkinds.Rd | 2 + man/mkinfit.Rd | 8 +-- man/mkinpredict.Rd | 108 ++++++++++++++-------------- test.log | 30 ++++---- tests/testthat/FOCUS_2006_D.csf | 2 +- vignettes/FOCUS_D.html | 120 ++++++++++++++++--------------- 12 files changed, 248 insertions(+), 226 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 99a9cc76..96f4637d 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -2,7 +2,7 @@ Package: mkin Type: Package Title: Kinetic Evaluation of Chemical Degradation Data Version: 0.9.49.11 -Date: 2020-04-20 +Date: 2020-05-20 Authors@R: c(person("Johannes", "Ranke", role = c("aut", "cre", "cph"), email = "jranke@uni-bremen.de", comment = c(ORCID = "0000-0003-4371-6538")), diff --git a/R/mkinfit.R b/R/mkinfit.R index 5c092612..2b7e71cb 100644 --- a/R/mkinfit.R +++ b/R/mkinfit.R @@ -185,13 +185,13 @@ if(getRversion() >= '2.15.1') utils::globalVariables(c("name", "time", "value")) #' # Fit the model to the FOCUS example dataset D using defaults #' print(system.time(fit <- mkinfit(SFO_SFO, FOCUS_2006_D, #' solution_type = "eigen", quiet = TRUE))) -#' coef(fit) +#' parms(fit) #' endpoints(fit) #' \dontrun{ #' # deSolve is slower when no C compiler (gcc) was available during model generation #' print(system.time(fit.deSolve <- mkinfit(SFO_SFO, FOCUS_2006_D, #' solution_type = "deSolve"))) -#' coef(fit.deSolve) +#' parms(fit.deSolve) #' endpoints(fit.deSolve) #' } #' @@ -926,6 +926,6 @@ mkinfit <- function(mkinmod, observed, fit$version <- as.character(utils::packageVersion("mkin")) fit$Rversion <- paste(R.version$major, R.version$minor, sep=".") - class(fit) <- c("mkinfit", "modFit") + class(fit) <- c("mkinfit") return(fit) } diff --git a/R/mkinmod.R b/R/mkinmod.R index ca1402fd..099e1155 100644 --- a/R/mkinmod.R +++ b/R/mkinmod.R @@ -421,6 +421,46 @@ mkinmod <- function(..., use_of_ff = "min", speclist = NULL, quiet = FALSE, verb } # }}} + # Attach a degradation function if an analytical solution is available {{{ + parent_type = spec[[1]]$type + parent_name = names(spec)[[1]] + if (length(spec) == 1) { + odeparm_map <- switch(parent_type, + SFO = c( + k = if(use_of_ff == "min") paste("k", parent_name, "sink", sep = "_") + else paste("k", parent_name, sep = "_")), + FOMC = c(alpha = "alpha", beta = "beta"), + IORE = c( + k__iore = if(use_of_ff == "min") paste("k__iore", parent_name, "sink", sep = "_") + else paste("k__iore", parent_name, sep = "_"), + N = paste("N", parent_name, sep = "_")), + DFOP = c(k1 = "k1", k2 = "k2", g = "g"), + HS = c(k1 = "k1", k2 = "k2", tb = "tb"), + SFORB = c( + k_12 = paste("k", parent_name, "free_bound", sep = "_"), + k_21 = paste("k", parent_name, "bound_free", sep = "_"), + k_1output = paste("k", parent_name, "free_sink", sep = "_")), + logistic = c(kmax = "kmax", k0 = "k0", r = "r") + ) + odeparm_rev_map <- names(odeparm_map) + names(odeparm_rev_map) <- odeparm_map + + model$deg_func <- function(odeini, odeparms, outtimes) { + parent_func <- getFromNamespace(paste0(parent_type, ".solution"), "mkin") + odeparm_list <- as.list(odeparms) + names(odeparm_list) <- odeparm_rev_map[names(odeparm_list)] + + values <- do.call(parent_func, + args = c( + list(t = outtimes, parent.0 = odeini[1]), + odeparm_list)) + out <- data.frame(outtimes, values) + names(out) <- c("time", parent_name) + return(out) + } + } + # }}} + class(model) <- "mkinmod" return(model) } diff --git a/R/mkinpredict.R b/R/mkinpredict.R index 16ee7903..0f8e83bb 100644 --- a/R/mkinpredict.R +++ b/R/mkinpredict.R @@ -41,57 +41,63 @@ #' @author Johannes Ranke #' @examples #' -#' SFO <- mkinmod(degradinol = mkinsub("SFO")) -#' # Compare solution types -#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, -#' solution_type = "analytical") -#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, -#' solution_type = "deSolve") -#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, -#' solution_type = "deSolve", use_compiled = FALSE) -#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, -#' solution_type = "eigen") +#' SFO <- mkinmod(degradinol = mkinsub("SFO")) +#' # Compare solution types +#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, +#' solution_type = "analytical") +#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, +#' solution_type = "deSolve") +#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, +#' solution_type = "deSolve", use_compiled = FALSE) +#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, +#' solution_type = "eigen") #' +#' # Compare integration methods to analytical solution +#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, +#' solution_type = "analytical")[21,] +#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, +#' method = "lsoda")[21,] +#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, +#' method = "ode45")[21,] +#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, +#' method = "rk4")[21,] +#' # rk4 is not as precise here #' -#' # Compare integration methods to analytical solution -#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, -#' solution_type = "analytical")[21,] -#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, -#' method = "lsoda")[21,] -#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, -#' method = "ode45")[21,] -#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, -#' method = "rk4")[21,] -#' # rk4 is not as precise here +#' # The number of output times used to make a lot of difference until the +#' # default for atol was adjusted +#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), +#' seq(0, 20, by = 0.1))[201,] +#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), +#' seq(0, 20, by = 0.01))[2001,] #' -#' # The number of output times used to make a lot of difference until the -#' # default for atol was adjusted -#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), -#' seq(0, 20, by = 0.1))[201,] -#' mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), -#' seq(0, 20, by = 0.01))[2001,] +#' # Check compiled model versions - they are faster than the eigenvalue based solutions! +#' SFO_SFO = mkinmod(parent = list(type = "SFO", to = "m1"), +#' m1 = list(type = "SFO")) +#' if(require(rbenchmark)) { +#' benchmark( +#' eigen = mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), +#' c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), +#' solution_type = "eigen")[201,], +#' deSolve_compiled = mkinpredict(SFO_SFO, +#' c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), +#' c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), +#' solution_type = "deSolve")[201,], +#' deSolve = mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), +#' c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), +#' solution_type = "deSolve", use_compiled = FALSE)[201,], +#' replications = 10) +#' } #' -#' # Check compiled model versions - they are faster than the eigenvalue based solutions! -#' SFO_SFO = mkinmod(parent = list(type = "SFO", to = "m1"), -#' m1 = list(type = "SFO")) -#' system.time( -#' print(mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), -#' c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), -#' solution_type = "eigen")[201,])) -#' system.time( -#' print(mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), -#' c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), -#' solution_type = "deSolve")[201,])) -#' system.time( -#' print(mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), -#' c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), -#' solution_type = "deSolve", use_compiled = FALSE)[201,])) +#' # Since mkin 0.9.49.11 we also have analytical solutions for some models, including SFO-SFO +#' # deSolve = mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), +#' # c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), +#' # solution_type = "analytical", use_compiled = FALSE)[201,], #' -#' \dontrun{ -#' # Predict from a fitted model -#' f <- mkinfit(SFO_SFO, FOCUS_2006_C) -#' head(mkinpredict(f)) -#' } +#' \dontrun{ +#' # Predict from a fitted model +#' f <- mkinfit(SFO_SFO, FOCUS_2006_C, quiet = TRUE) +#' head(mkinpredict(f)) +#' } #' #' @export mkinpredict <- function(x, odeparms, odeini, @@ -124,54 +130,17 @@ mkinpredict.mkinmod <- function(x, odeini <- odeini[mod_vars] } - # Create function for evaluation of expressions with ode parameters and initial values - evalparse <- function(string) - { - eval(parse(text=string), as.list(c(odeparms, odeini))) - } - - # Create a function calculating the differentials specified by the model - # if necessary if (solution_type == "analytical") { - parent.type = names(x$map[[1]])[1] - parent.name = names(x$diffs)[[1]] - o <- switch(parent.type, - SFO = SFO.solution(outtimes, - evalparse(parent.name), - ifelse(x$use_of_ff == "min", - evalparse(paste("k", parent.name, "sink", sep="_")), - evalparse(paste("k", parent.name, sep="_")))), - FOMC = FOMC.solution(outtimes, - evalparse(parent.name), - evalparse("alpha"), evalparse("beta")), - IORE = IORE.solution(outtimes, - evalparse(parent.name), - ifelse(x$use_of_ff == "min", - evalparse(paste("k__iore", parent.name, "sink", sep="_")), - evalparse(paste("k__iore", parent.name, sep="_"))), - evalparse("N_parent")), - DFOP = DFOP.solution(outtimes, - evalparse(parent.name), - evalparse("k1"), evalparse("k2"), - evalparse("g")), - HS = HS.solution(outtimes, - evalparse(parent.name), - evalparse("k1"), evalparse("k2"), - evalparse("tb")), - SFORB = SFORB.solution(outtimes, - evalparse(parent.name), - evalparse(paste("k", parent.name, "bound", sep="_")), - evalparse(paste("k", sub("free", "bound", parent.name), "free", sep="_")), - evalparse(paste("k", parent.name, "sink", sep="_"))), - logistic = logistic.solution(outtimes, - evalparse(parent.name), - evalparse("kmax"), evalparse("k0"), - evalparse("r")) - ) - out <- data.frame(outtimes, o) - names(out) <- c("time", sub("_free", "", parent.name)) + out <- x$deg_func(odeini = odeini, + odeparms = odeparms, outtimes = outtimes) } + if (solution_type == "eigen") { + + evalparse <- function(string) { + eval(parse(text=string), as.list(c(odeparms, odeini))) + } + coefmat.num <- matrix(sapply(as.vector(x$coefmat), evalparse), nrow = length(mod_vars)) e <- eigen(coefmat.num) @@ -184,6 +153,7 @@ mkinpredict.mkinmod <- function(x, out <- data.frame(outtimes, t(o)) names(out) <- c("time", mod_vars) } + if (solution_type == "deSolve") { if (!is.null(x$cf) & use_compiled[1] != FALSE) { out <- ode( diff --git a/build.log b/build.log index 4fd091b2..e1e92209 100644 --- a/build.log +++ b/build.log @@ -5,6 +5,5 @@ * creating vignettes ... OK * checking for LF line-endings in source and make files and shell scripts * checking for empty or unneeded directories -* looking to see if a ‘data/datalist’ file should be added * building ‘mkin_0.9.49.11.tar.gz’ diff --git a/check.log b/check.log index df260cca..ddf3c6f4 100644 --- a/check.log +++ b/check.log @@ -1,5 +1,5 @@ * using log directory ‘/home/jranke/git/mkin/mkin.Rcheck’ -* using R version 3.6.3 (2020-02-29) +* using R version 4.0.0 (2020-04-24) * using platform: x86_64-pc-linux-gnu (64-bit) * using session charset: UTF-8 * using options ‘--no-tests --as-cran’ @@ -63,6 +63,7 @@ Maintainer: ‘Johannes Ranke ’ * checking package vignettes in ‘inst/doc’ ... OK * checking re-building of vignette outputs ... OK * checking PDF version of manual ... OK +* checking for non-standard things in the check directory ... OK * checking for detritus in the temp directory ... OK * DONE diff --git a/man/mkinds.Rd b/man/mkinds.Rd index 3bbb1c4b..8d1860e8 100644 --- a/man/mkinds.Rd +++ b/man/mkinds.Rd @@ -47,6 +47,7 @@ and value in order to be compatible with mkinfit} } \if{html}{\out{
}} \if{html}{\out{}} +\if{latex}{\out{\hypertarget{method-new}{}}} \subsection{Method \code{new()}}{ Create a new mkinds object \subsection{Usage}{ @@ -69,6 +70,7 @@ Create a new mkinds object } \if{html}{\out{
}} \if{html}{\out{}} +\if{latex}{\out{\hypertarget{method-clone}{}}} \subsection{Method \code{clone()}}{ The objects of this class are cloneable with this method. \subsection{Usage}{ diff --git a/man/mkinfit.Rd b/man/mkinfit.Rd index 45036361..89d2ddd6 100644 --- a/man/mkinfit.Rd +++ b/man/mkinfit.Rd @@ -28,8 +28,8 @@ mkinfit( rtol = 1e-10, n.outtimes = 100, error_model = c("const", "obs", "tc"), - error_model_algorithm = c("auto", "d_3", "direct", "twostep", "threestep", "fourstep", - "IRLS", "OLS"), + error_model_algorithm = c("auto", "d_3", "direct", "twostep", "threestep", + "fourstep", "IRLS", "OLS"), reweight.tol = 1e-08, reweight.max.iter = 10, trace_parms = FALSE, @@ -238,13 +238,13 @@ SFO_SFO <- mkinmod( # Fit the model to the FOCUS example dataset D using defaults print(system.time(fit <- mkinfit(SFO_SFO, FOCUS_2006_D, solution_type = "eigen", quiet = TRUE))) -coef(fit) +parms(fit) endpoints(fit) \dontrun{ # deSolve is slower when no C compiler (gcc) was available during model generation print(system.time(fit.deSolve <- mkinfit(SFO_SFO, FOCUS_2006_D, solution_type = "deSolve"))) -coef(fit.deSolve) +parms(fit.deSolve) endpoints(fit.deSolve) } diff --git a/man/mkinpredict.Rd b/man/mkinpredict.Rd index 17d7ef21..366d5b83 100644 --- a/man/mkinpredict.Rd +++ b/man/mkinpredict.Rd @@ -100,57 +100,63 @@ kinetic parameters and initial values for the state variables. } \examples{ - SFO <- mkinmod(degradinol = mkinsub("SFO")) - # Compare solution types - mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, - solution_type = "analytical") - mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, - solution_type = "deSolve") - mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, - solution_type = "deSolve", use_compiled = FALSE) - mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, - solution_type = "eigen") - - - # Compare integration methods to analytical solution - mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, - solution_type = "analytical")[21,] - mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, - method = "lsoda")[21,] - mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, - method = "ode45")[21,] - mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, - method = "rk4")[21,] - # rk4 is not as precise here - - # The number of output times used to make a lot of difference until the - # default for atol was adjusted - mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), - seq(0, 20, by = 0.1))[201,] - mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), - seq(0, 20, by = 0.01))[2001,] - - # Check compiled model versions - they are faster than the eigenvalue based solutions! - SFO_SFO = mkinmod(parent = list(type = "SFO", to = "m1"), - m1 = list(type = "SFO")) - system.time( - print(mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), - c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), - solution_type = "eigen")[201,])) - system.time( - print(mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), - c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), - solution_type = "deSolve")[201,])) - system.time( - print(mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), - c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), - solution_type = "deSolve", use_compiled = FALSE)[201,])) - - \dontrun{ - # Predict from a fitted model - f <- mkinfit(SFO_SFO, FOCUS_2006_C) - head(mkinpredict(f)) - } +SFO <- mkinmod(degradinol = mkinsub("SFO")) +# Compare solution types +mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, + solution_type = "analytical") +mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, + solution_type = "deSolve") +mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, + solution_type = "deSolve", use_compiled = FALSE) +mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, + solution_type = "eigen") + +# Compare integration methods to analytical solution +mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, + solution_type = "analytical")[21,] +mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, + method = "lsoda")[21,] +mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, + method = "ode45")[21,] +mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20, + method = "rk4")[21,] +# rk4 is not as precise here + +# The number of output times used to make a lot of difference until the +# default for atol was adjusted +mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), + seq(0, 20, by = 0.1))[201,] +mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), + seq(0, 20, by = 0.01))[2001,] + +# Check compiled model versions - they are faster than the eigenvalue based solutions! +SFO_SFO = mkinmod(parent = list(type = "SFO", to = "m1"), + m1 = list(type = "SFO")) +if(require(rbenchmark)) { + benchmark( + eigen = mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), + c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), + solution_type = "eigen")[201,], + deSolve_compiled = mkinpredict(SFO_SFO, + c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), + c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), + solution_type = "deSolve")[201,], + deSolve = mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), + c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), + solution_type = "deSolve", use_compiled = FALSE)[201,], + replications = 10) +} + +# Since mkin 0.9.49.11 we also have analytical solutions for some models, including SFO-SFO +# deSolve = mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01), +# c(parent = 100, m1 = 0), seq(0, 20, by = 0.1), +# solution_type = "analytical", use_compiled = FALSE)[201,], + +\dontrun{ + # Predict from a fitted model + f <- mkinfit(SFO_SFO, FOCUS_2006_C, quiet = TRUE) + head(mkinpredict(f)) +} } \author{ diff --git a/test.log b/test.log index a7879c33..b57bd69a 100644 --- a/test.log +++ b/test.log @@ -2,32 +2,32 @@ Loading mkin Testing mkin ✔ | OK F W S | Context ✔ | 2 | Export dataset for reading into CAKE -✔ | 13 | Results for FOCUS D established in expertise for UBA (Ranke 2014) [3.5 s] -✔ | 4 | Calculation of FOCUS chi2 error levels [2.1 s] -✔ | 6 | Fitting the SFORB model [8.8 s] +✔ | 13 | Results for FOCUS D established in expertise for UBA (Ranke 2014) [3.6 s] +✔ | 4 | Calculation of FOCUS chi2 error levels [2.2 s] +✔ | 6 | Fitting the SFORB model [9.1 s] ✔ | 5 | Calculation of Akaike weights -✔ | 10 | Confidence intervals and p-values [9.4 s] -✔ | 14 | Error model fitting [37.2 s] +✔ | 10 | Confidence intervals and p-values [9.6 s] +✔ | 14 | Error model fitting [38.5 s] ✔ | 6 | Test fitting the decline of metabolites from their maximum [0.8 s] -✔ | 1 | Fitting the logistic model [0.9 s] +✔ | 1 | Fitting the logistic model [0.8 s] ✔ | 1 | Test dataset class mkinds used in gmkin ✔ | 12 | Special cases of mkinfit calls [2.3 s] ✔ | 8 | mkinmod model generation and printing [0.2 s] -✔ | 3 | Model predictions with mkinpredict [0.3 s] -✔ | 16 | Evaluations according to 2015 NAFTA guidance [4.1 s] -✔ | 9 | Nonlinear mixed-effects models [11.9 s] -✔ | 4 | Calculation of maximum time weighted average concentrations (TWAs) [2.2 s] +✔ | 3 | Model predictions with mkinpredict [0.4 s] +✔ | 16 | Evaluations according to 2015 NAFTA guidance [4.2 s] +✔ | 9 | Nonlinear mixed-effects models [12.3 s] +✔ | 4 | Calculation of maximum time weighted average concentrations (TWAs) [2.4 s] ✔ | 3 | Summary -✔ | 14 | Plotting [4.9 s] +✔ | 14 | Plotting [5.1 s] ✔ | 4 | AIC calculation ✔ | 4 | Residuals extracted from mkinfit models -✔ | 2 | Complex test case from Schaefer et al. (2007) Piacenza paper [5.4 s] +✔ | 2 | Complex test case from Schaefer et al. (2007) Piacenza paper [5.3 s] ✔ | 1 | Summaries of old mkinfit objects -✔ | 4 | Results for synthetic data established in expertise for UBA (Ranke 2014) [6.9 s] -✔ | 9 | Hypothesis tests [36.7 s] +✔ | 4 | Results for synthetic data established in expertise for UBA (Ranke 2014) [7.2 s] +✔ | 9 | Hypothesis tests [38.2 s] ══ Results ═════════════════════════════════════════════════════════════════════ -Duration: 137.8 s +Duration: 142.5 s OK: 155 Failed: 0 diff --git a/tests/testthat/FOCUS_2006_D.csf b/tests/testthat/FOCUS_2006_D.csf index 8faa19db..fa96de74 100644 --- a/tests/testthat/FOCUS_2006_D.csf +++ b/tests/testthat/FOCUS_2006_D.csf @@ -5,7 +5,7 @@ Description: MeasurementUnits: % AR TimeUnits: days Comments: Created using mkin::CAKE_export -Date: 2020-04-22 +Date: 2020-05-06 Optimiser: IRLS [Data] diff --git a/vignettes/FOCUS_D.html b/vignettes/FOCUS_D.html index 7a12d221..ff9d4596 100644 --- a/vignettes/FOCUS_D.html +++ b/vignettes/FOCUS_D.html @@ -1,17 +1,17 @@ - + - + - + Example evaluation of FOCUS Example Dataset D @@ -270,7 +270,6 @@ code { } img { max-width:100%; - height: auto; } .tabbed-pane { padding-top: 12px; @@ -332,6 +331,7 @@ summary { border: none; display: inline-block; border-radius: 4px; + background-color: transparent; } .tabset-dropdown > .nav-tabs.nav-tabs-open > li { @@ -344,18 +344,6 @@ summary { } - - @@ -377,7 +365,7 @@ $(document).ready(function () {

Example evaluation of FOCUS Example Dataset D

Johannes Ranke

-

2019-05-02

+

2020-05-06

@@ -441,20 +429,20 @@ print(FOCUS_2006_D) ## "d_m1 = + k_parent_m1 * parent - k_m1_sink * m1"

We do the fitting without progress report (quiet = TRUE).

fit <- mkinfit(SFO_SFO, FOCUS_2006_D, quiet = TRUE)
-
## Warning in mkinfit(SFO_SFO, FOCUS_2006_D, quiet = TRUE): Observations with
-## value of zero were removed from the data
+
## Warning in mkinfit(SFO_SFO, FOCUS_2006_D, quiet = TRUE): Observations with value
+## of zero were removed from the data

A plot of the fit including a residual plot for both observed variables is obtained using the plot_sep method for mkinfit objects, which shows separate graphs for all compounds and their residuals.

plot_sep(fit, lpos = c("topright", "bottomright"))
-

+

Confidence intervals for the parameter estimates are obtained using the mkinparplot function.

mkinparplot(fit)
-

+

A comprehensive report of the results is obtained using the summary method for mkinfit objects.

summary(fit)
-
## mkin version used for fitting:    0.9.49.4 
-## R version used for fitting:       3.6.0 
-## Date of fit:     Thu May  2 18:43:48 2019 
-## Date of summary: Thu May  2 18:43:49 2019 
+
## mkin version used for fitting:    0.9.49.11 
+## R version used for fitting:       4.0.0 
+## Date of fit:     Wed May  6 19:54:54 2020 
+## Date of summary: Wed May  6 19:54:54 2020 
 ## 
 ## Equations:
 ## d_parent/dt = - k_parent_sink * parent - k_parent_m1 * parent
@@ -462,18 +450,18 @@ print(FOCUS_2006_D)

 ## 
 ## Model predictions using solution type deSolve 
 ## 
-## Fitted using 396 model solutions performed in 1.002 s
+## Fitted using 389 model solutions performed in 1.041 s
 ## 
-## Error model:
-## Constant variance 
+## Error model: Constant variance 
+## 
+## Error model algorithm: OLS 
 ## 
 ## Starting values for parameters to be optimised:
-##                    value   type
-## parent_0      100.750000  state
-## k_parent_sink   0.100000 deparm
-## k_parent_m1     0.100100 deparm
-## k_m1_sink       0.100200 deparm
-## sigma           3.125504  error
+##                  value   type
+## parent_0      100.7500  state
+## k_parent_sink   0.1000 deparm
+## k_parent_m1     0.1001 deparm
+## k_m1_sink       0.1002 deparm
 ## 
 ## Starting values for the transformed parameters actually optimised:
 ##                        value lower upper
@@ -481,7 +469,6 @@ print(FOCUS_2006_D)
## log_k_parent_sink -2.302585 -Inf Inf ## log_k_parent_m1 -2.301586 -Inf Inf ## log_k_m1_sink -2.300587 -Inf Inf -## sigma 3.125504 0 Inf ## ## Fixed parameter values: ## value type @@ -496,18 +483,18 @@ print(FOCUS_2006_D) ## sigma 3.126 0.35850 2.396 3.855 ## ## Parameter correlation: -## parent_0 log_k_parent_sink log_k_parent_m1 -## parent_0 1.000e+00 6.067e-01 -6.372e-02 -## log_k_parent_sink 6.067e-01 1.000e+00 -8.550e-02 -## log_k_parent_m1 -6.372e-02 -8.550e-02 1.000e+00 -## log_k_m1_sink -1.688e-01 -6.252e-01 4.731e-01 -## sigma -3.368e-07 3.365e-08 8.420e-07 -## log_k_m1_sink sigma -## parent_0 -1.688e-01 -3.368e-07 -## log_k_parent_sink -6.252e-01 3.365e-08 -## log_k_parent_m1 4.731e-01 8.420e-07 -## log_k_m1_sink 1.000e+00 1.958e-08 -## sigma 1.958e-08 1.000e+00 +## parent_0 log_k_parent_sink log_k_parent_m1 log_k_m1_sink +## parent_0 1.000e+00 6.067e-01 -6.372e-02 -1.688e-01 +## log_k_parent_sink 6.067e-01 1.000e+00 -8.550e-02 -6.252e-01 +## log_k_parent_m1 -6.372e-02 -8.550e-02 1.000e+00 4.731e-01 +## log_k_m1_sink -1.688e-01 -6.252e-01 4.731e-01 1.000e+00 +## sigma 5.287e-10 3.306e-09 4.421e-08 -3.319e-10 +## sigma +## parent_0 5.287e-10 +## log_k_parent_sink 3.306e-09 +## log_k_parent_m1 4.421e-08 +## log_k_m1_sink -3.319e-10 +## sigma 1.000e+00 ## ## Backtransformed parameters: ## Confidence intervals for internally transformed parameters are asymmetric. @@ -543,10 +530,10 @@ print(FOCUS_2006_D) ## 0 parent 102.04 99.59848 2.442e+00 ## 1 parent 93.50 90.23787 3.262e+00 ## 1 parent 92.50 90.23787 2.262e+00 -## 3 parent 63.23 74.07319 -1.084e+01 -## 3 parent 68.99 74.07319 -5.083e+00 -## 7 parent 52.32 49.91206 2.408e+00 -## 7 parent 55.13 49.91206 5.218e+00 +## 3 parent 63.23 74.07320 -1.084e+01 +## 3 parent 68.99 74.07320 -5.083e+00 +## 7 parent 52.32 49.91207 2.408e+00 +## 7 parent 55.13 49.91207 5.218e+00 ## 14 parent 27.27 25.01257 2.257e+00 ## 14 parent 26.64 25.01257 1.627e+00 ## 21 parent 11.50 12.53462 -1.035e+00 @@ -556,23 +543,23 @@ print(FOCUS_2006_D) ## 50 parent 0.69 0.71624 -2.624e-02 ## 50 parent 0.63 0.71624 -8.624e-02 ## 75 parent 0.05 0.06074 -1.074e-02 -## 75 parent 0.06 0.06074 -7.381e-04 +## 75 parent 0.06 0.06074 -7.382e-04 ## 1 m1 4.84 4.80296 3.704e-02 ## 1 m1 5.64 4.80296 8.370e-01 ## 3 m1 12.91 13.02400 -1.140e-01 ## 3 m1 12.96 13.02400 -6.400e-02 ## 7 m1 22.97 25.04476 -2.075e+00 ## 7 m1 24.47 25.04476 -5.748e-01 -## 14 m1 41.69 36.69003 5.000e+00 -## 14 m1 33.21 36.69003 -3.480e+00 +## 14 m1 41.69 36.69002 5.000e+00 +## 14 m1 33.21 36.69002 -3.480e+00 ## 21 m1 44.37 41.65310 2.717e+00 ## 21 m1 46.44 41.65310 4.787e+00 ## 35 m1 41.22 43.31312 -2.093e+00 ## 35 m1 37.95 43.31312 -5.363e+00 -## 50 m1 41.19 41.21832 -2.832e-02 -## 50 m1 40.01 41.21832 -1.208e+00 -## 75 m1 40.09 36.44704 3.643e+00 -## 75 m1 33.85 36.44704 -2.597e+00 +## 50 m1 41.19 41.21831 -2.831e-02 +## 50 m1 40.01 41.21831 -1.208e+00 +## 75 m1 40.09 36.44703 3.643e+00 +## 75 m1 33.85 36.44703 -2.597e+00 ## 100 m1 31.04 31.98163 -9.416e-01 ## 100 m1 33.13 31.98163 1.148e+00 ## 120 m1 25.15 28.78984 -3.640e+00 @@ -596,6 +583,23 @@ $(document).ready(function () { + + + + + + +