diff options
| -rw-r--r-- | NAMESPACE | 3 | ||||
| -rw-r--r-- | NEWS.md | 4 | ||||
| -rw-r--r-- | R/llhist.R | 3 | ||||
| -rw-r--r-- | R/multistart.R | 6 | ||||
| -rw-r--r-- | R/parplot.R (renamed from R/parhist.R) | 12 | ||||
| -rw-r--r-- | R/saem.R | 2 | ||||
| -rw-r--r-- | log/test.log | 42 | ||||
| -rw-r--r-- | man/logLik.saem.mmkin.Rd | 2 | ||||
| -rw-r--r-- | man/multistart.Rd | 6 | ||||
| -rw-r--r-- | man/parplot.Rd (renamed from man/parhist.Rd) | 17 | ||||
| -rw-r--r-- | tests/testthat/_snaps/multistart/llhist-for-biphasic-saemix-fit.svg | 1 | ||||
| -rw-r--r-- | tests/testthat/_snaps/multistart/llhist-for-sfo-fit.svg | 2 | ||||
| -rw-r--r-- | tests/testthat/_snaps/multistart/parplot-for-biphasic-saemix-fit.svg (renamed from tests/testthat/_snaps/multistart/parhist-for-biphasic-saemix-fit.svg) | 0 | ||||
| -rw-r--r-- | tests/testthat/_snaps/multistart/parplot-for-sfo-fit.svg (renamed from tests/testthat/_snaps/multistart/parhist-for-sfo-fit.svg) | 0 | ||||
| -rw-r--r-- | tests/testthat/test_multistart.R | 8 | ||||
| -rw-r--r-- | vignettes/web_only/multistart.rmd | 6 |
16 files changed, 63 insertions, 51 deletions
@@ -36,6 +36,7 @@ S3method(parms,mkinfit) S3method(parms,mmkin) S3method(parms,multistart) S3method(parms,saem.mmkin) +S3method(parplot,multistart.saem.mmkin) S3method(plot,mixed.mmkin) S3method(plot,mkinfit) S3method(plot,mmkin) @@ -126,8 +127,8 @@ export(nafta) export(nlme) export(nlme_data) export(nlme_function) -export(parhist) export(parms) +export(parplot) export(plot_err) export(plot_res) export(plot_sep) @@ -1,6 +1,6 @@ # mkin 1.1.2 -- 'R/multistart.R': New method for testing multiple start parameters for hierarchical model fits, with diagnostic plotting functions 'llhist' and 'parhist'. +- 'R/multistart.R': New method for testing multiple start parameters for hierarchical model fits, with diagnostic plotting functions 'llhist' and 'parplot'. - 'R/mhmkin.R': New method for performing multiple hierarchical mkin fits in one function call, optionally in parallel. @@ -8,7 +8,7 @@ - 'R/saem.R': 'logLik' and 'update' methods for 'saem.mmkin' objects. -- 'R/convergence.R': New generic to show convergence information with methods for 'mmkin' and 'mhmkin' objects. +- 'R/status.R': New generic to show status information for fit array objects with methods for 'mmkin' and 'mhmkin' objects. - 'R/illparms.R': New generic to show ill-defined parameters with methods for 'mkinfit', 'mmkin', 'saem.mmkin' and 'mhmkin' objects. @@ -25,6 +25,7 @@ llhist <- function(object, breaks = "Sturges", lpos = "topleft", main = "", stop("llhist is only implemented for multistart.saem.mmkin objects") } + ll_orig <- logLik(attr(object, "orig")) ll <- stats::na.omit(sapply(object, llfunc)) par(las = 1) @@ -34,7 +35,7 @@ llhist <- function(object, breaks = "Sturges", lpos = "topleft", main = "", freq_factor <- h$counts[1] / h$density[1] - abline(v = logLik(attr(object, "orig")), col = 2) + abline(v = ll_orig, col = 2) legend(lpos, inset = c(0.05, 0.05), bty = "n", lty = 1, col = c(2), diff --git a/R/multistart.R b/R/multistart.R index 29ccdc44..14683c11 100644 --- a/R/multistart.R +++ b/R/multistart.R @@ -17,7 +17,7 @@ #' @param x The multistart object to print #' @return A list of [saem.mmkin] objects, with class attributes #' 'multistart.saem.mmkin' and 'multistart'. -#' @seealso [parhist], [llhist] +#' @seealso [parplot], [llhist] #' #' @references Duchesne R, Guillemin A, Gandrillon O, Crauste F. Practical #' identifiability in the frame of nonlinear mixed effects models: the example @@ -40,7 +40,7 @@ #' f_mmkin <- mmkin("DFOP", dmta_ds, error_model = "tc", cores = 7, quiet = TRUE) #' f_saem_full <- saem(f_mmkin) #' f_saem_full_multi <- multistart(f_saem_full, n = 16, cores = 16) -#' parhist(f_saem_full_multi, lpos = "bottomleft") +#' parplot(f_saem_full_multi, lpos = "bottomleft") #' illparms(f_saem_full) #' #' f_saem_reduced <- update(f_saem_full, no_random_effect = "log_k2") @@ -52,7 +52,7 @@ #' cl <- makePSOCKcluster(12) #' clusterExport(cl, "f_mmkin") #' f_saem_reduced_multi <- multistart(f_saem_reduced, n = 16, cluster = cl) -#' parhist(f_saem_reduced_multi, lpos = "topright") +#' parplot(f_saem_reduced_multi, lpos = "topright") #' } multistart <- function(object, n = 50, cores = if (Sys.info()["sysname"] == "Windows") 1 else parallel::detectCores(), diff --git a/R/parhist.R b/R/parplot.R index 0f9c9964..627a4eb8 100644 --- a/R/parhist.R +++ b/R/parplot.R @@ -1,4 +1,4 @@ -#' Plot parameter distributions from multistart objects +#' Plot parameter variability of multistart objects #' #' Produces a boxplot with all parameters from the multiple runs, scaled #' either by the parameters of the run with the highest likelihood, @@ -18,7 +18,13 @@ #' @seealso [multistart] #' @importFrom stats median #' @export -parhist <- function(object, llmin = -Inf, scale = c("best", "median"), +parplot <- function(object, ...) { + UseMethod("parplot") +} + +#' @rdname parplot +#' @export +parplot.multistart.saem.mmkin <- function(object, llmin = -Inf, scale = c("best", "median"), lpos = "bottomleft", main = "", ...) { oldpar <- par(no.readonly = TRUE) @@ -35,7 +41,7 @@ parhist <- function(object, llmin = -Inf, scale = c("best", "median"), else return(logLik(object$so)) } } else { - stop("parhist is only implemented for multistart.saem.mmkin objects") + stop("parplot is only implemented for multistart.saem.mmkin objects") } ll <- sapply(object, llfunc) selected <- which(ll > llmin) @@ -726,7 +726,7 @@ saemix_data <- function(object, covariates = NULL, verbose = FALSE, ...) { #' @param \dots Passed to [saemix::logLik.SaemixObject] #' @param method Passed to [saemix::logLik.SaemixObject] #' @export -logLik.saem.mmkin <- function(object, ..., method = c("lin", "is", "gq")) { +logLik.saem.mmkin <- function(object, ..., method = c("is", "lin", "gq")) { method <- match.arg(method) return(logLik(object$so, method = method)) } diff --git a/log/test.log b/log/test.log index 429a2e02..0e2ca6b2 100644 --- a/log/test.log +++ b/log/test.log @@ -1,53 +1,53 @@ ℹ Testing mkin ✔ | F W S OK | Context ✔ | 5 | AIC calculation -✔ | 5 | Analytical solutions for coupled models [3.3s] +✔ | 5 | Analytical solutions for coupled models [3.4s] ✔ | 5 | Calculation of Akaike weights ✔ | 3 | Export dataset for reading into CAKE ✔ | 12 | Confidence intervals and p-values [1.0s] -✔ | 1 12 | Dimethenamid data from 2018 [43.2s] +✔ | 1 12 | Dimethenamid data from 2018 [63.1s] ──────────────────────────────────────────────────────────────────────────────── Skip (test_dmta.R:98: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 [11.5s] +✔ | 14 | Error model fitting [5.5s] ✔ | 5 | Time step normalisation -✔ | 4 | Calculation of FOCUS chi2 error levels [1.3s] -✔ | 14 | Results for FOCUS D established in expertise for UBA (Ranke 2014) [1.9s] -✔ | 4 | Test fitting the decline of metabolites from their maximum [0.9s] -✔ | 1 | Fitting the logistic model [0.6s] -✔ | 7 | Batch fitting and diagnosing hierarchical kinetic models [30.9s] -✔ | 1 12 | Nonlinear mixed-effects models [0.4s] +✔ | 4 | Calculation of FOCUS chi2 error levels [0.6s] +✔ | 14 | Results for FOCUS D established in expertise for UBA (Ranke 2014) [0.8s] +✔ | 4 | Test fitting the decline of metabolites from their maximum [0.4s] +✔ | 1 | Fitting the logistic model [0.2s] +✔ | 7 | Batch fitting and diagnosing hierarchical kinetic models [28.6s] +✔ | 1 12 | Nonlinear mixed-effects models [0.6s] ──────────────────────────────────────────────────────────────────────────────── Skip (test_mixed.R:74:3): saemix results are reproducible for biphasic fits Reason: Fitting with saemix takes around 10 minutes when using deSolve ──────────────────────────────────────────────────────────────────────────────── ✔ | 3 | Test dataset classes mkinds and mkindsg -✔ | 10 | Special cases of mkinfit calls [0.7s] -✔ | 3 | mkinfit features [1.0s] -✔ | 8 | mkinmod model generation and printing [0.3s] -✔ | 3 | Model predictions with mkinpredict [0.3s] -✔ | 7 | Multistart method for saem.mmkin models [39.8s] -✔ | 16 | Evaluations according to 2015 NAFTA guidance [2.6s] -✔ | 9 | Nonlinear mixed-effects models with nlme [8.8s] -✔ | 16 | Plotting [10.1s] +✔ | 10 | Special cases of mkinfit calls [1.1s] +✔ | 3 | mkinfit features [1.7s] +✔ | 8 | mkinmod model generation and printing [0.4s] +✔ | 3 | Model predictions with mkinpredict [0.6s] +✔ | 7 | Multistart method for saem.mmkin models [93.9s] +✔ | 16 | Evaluations according to 2015 NAFTA guidance [7.0s] +✔ | 9 | Nonlinear mixed-effects models with nlme [13.0s] +✔ | 16 | Plotting [11.3s] ✔ | 4 | Residuals extracted from mkinfit models -✔ | 1 36 | saemix parent models [72.3s] +✔ | 1 36 | saemix parent models [72.6s] ──────────────────────────────────────────────────────────────────────────────── Skip (test_saemix_parent.R:152:3): We can also use mkin solution methods for saem Reason: This still takes almost 2.5 minutes although we do not solve ODEs ──────────────────────────────────────────────────────────────────────────────── ✔ | 2 | Complex test case from Schaefer et al. (2007) Piacenza paper [1.4s] ✔ | 11 | Processing of residue series -✔ | 7 | Fitting the SFORB model [3.8s] +✔ | 7 | Fitting the SFORB model [3.7s] ✔ | 1 | Summaries of old mkinfit objects ✔ | 5 | Summary [0.2s] ✔ | 4 | Results for synthetic data established in expertise for UBA (Ranke 2014) [2.2s] -✔ | 9 | Hypothesis tests [8.0s] +✔ | 9 | Hypothesis tests [8.1s] ✔ | 4 | Calculation of maximum time weighted average concentrations (TWAs) [2.2s] ══ Results ═════════════════════════════════════════════════════════════════════ -Duration: 249.3 s +Duration: 324.3 s ── Skipped tests ────────────────────────────────────────────────────────────── • Fitting this ODE model with saemix takes about 15 minutes on my system (1) diff --git a/man/logLik.saem.mmkin.Rd b/man/logLik.saem.mmkin.Rd index 603f4607..bd0bb72e 100644 --- a/man/logLik.saem.mmkin.Rd +++ b/man/logLik.saem.mmkin.Rd @@ -4,7 +4,7 @@ \alias{logLik.saem.mmkin} \title{logLik method for saem.mmkin objects} \usage{ -\method{logLik}{saem.mmkin}(object, ..., method = c("lin", "is", "gq")) +\method{logLik}{saem.mmkin}(object, ..., method = c("is", "lin", "gq")) } \arguments{ \item{object}{The fitted \link{saem.mmkin} object} diff --git a/man/multistart.Rd b/man/multistart.Rd index ebcc7d80..1f6773fe 100644 --- a/man/multistart.Rd +++ b/man/multistart.Rd @@ -77,7 +77,7 @@ dmta_ds[["Elliot 1"]] <- dmta_ds[["Elliot 2"]] <- NULL f_mmkin <- mmkin("DFOP", dmta_ds, error_model = "tc", cores = 7, quiet = TRUE) f_saem_full <- saem(f_mmkin) f_saem_full_multi <- multistart(f_saem_full, n = 16, cores = 16) -parhist(f_saem_full_multi, lpos = "bottomleft") +parplot(f_saem_full_multi, lpos = "bottomleft") illparms(f_saem_full) f_saem_reduced <- update(f_saem_full, no_random_effect = "log_k2") @@ -89,7 +89,7 @@ library(parallel) cl <- makePSOCKcluster(12) clusterExport(cl, "f_mmkin") f_saem_reduced_multi <- multistart(f_saem_reduced, n = 16, cluster = cl) -parhist(f_saem_reduced_multi, lpos = "topright") +parplot(f_saem_reduced_multi, lpos = "topright") } } \references{ @@ -99,5 +99,5 @@ of the in vitro erythropoiesis. BMC Bioinformatics. 2021 Oct 4;22(1):478. doi: 10.1186/s12859-021-04373-4. } \seealso{ -\link{parhist}, \link{llhist} +\link{parplot}, \link{llhist} } diff --git a/man/parhist.Rd b/man/parplot.Rd index f86ff734..37c5841d 100644 --- a/man/parhist.Rd +++ b/man/parplot.Rd @@ -1,10 +1,13 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/parhist.R -\name{parhist} -\alias{parhist} -\title{Plot parameter distributions from multistart objects} +% Please edit documentation in R/parplot.R +\name{parplot} +\alias{parplot} +\alias{parplot.multistart.saem.mmkin} +\title{Plot parameter variability of multistart objects} \usage{ -parhist( +parplot(object, ...) + +\method{parplot}{multistart.saem.mmkin}( object, llmin = -Inf, scale = c("best", "median"), @@ -16,6 +19,8 @@ parhist( \arguments{ \item{object}{The \link{multistart} object} +\item{\dots}{Passed to \link{boxplot}} + \item{llmin}{The minimum likelihood of objects to be shown} \item{scale}{By default, scale parameters using the best available fit. @@ -24,8 +29,6 @@ If 'median', parameters are scaled using the median parameters from all fits.} \item{lpos}{Positioning of the legend.} \item{main}{Title of the plot} - -\item{\dots}{Passed to \link{boxplot}} } \description{ Produces a boxplot with all parameters from the multiple runs, scaled diff --git a/tests/testthat/_snaps/multistart/llhist-for-biphasic-saemix-fit.svg b/tests/testthat/_snaps/multistart/llhist-for-biphasic-saemix-fit.svg index fa92c92a..fe38865d 100644 --- a/tests/testthat/_snaps/multistart/llhist-for-biphasic-saemix-fit.svg +++ b/tests/testthat/_snaps/multistart/llhist-for-biphasic-saemix-fit.svg @@ -50,6 +50,7 @@ <rect x='228.40' y='75.47' width='146.00' height='410.67' style='stroke-width: 0.75; fill: #D3D3D3;' /> <rect x='374.40' y='212.36' width='146.00' height='273.78' style='stroke-width: 0.75; fill: #D3D3D3;' /> <rect x='520.40' y='212.36' width='146.00' height='273.78' style='stroke-width: 0.75; fill: #D3D3D3;' /> +<line x1='352.00' y1='502.56' x2='352.00' y2='59.04' style='stroke-width: 0.75; stroke: #DF536B;' /> <line x1='101.38' y1='95.62' x2='122.98' y2='95.62' style='stroke-width: 0.75; stroke: #DF536B;' /> <text x='133.78' y='99.75' style='font-size: 12.00px; font-family: sans;' textLength='51.35px' lengthAdjust='spacingAndGlyphs'>original fit</text> </g> diff --git a/tests/testthat/_snaps/multistart/llhist-for-sfo-fit.svg b/tests/testthat/_snaps/multistart/llhist-for-sfo-fit.svg index 4342df9b..98513d06 100644 --- a/tests/testthat/_snaps/multistart/llhist-for-sfo-fit.svg +++ b/tests/testthat/_snaps/multistart/llhist-for-sfo-fit.svg @@ -50,7 +50,7 @@ <rect x='228.40' y='75.47' width='146.00' height='410.67' style='stroke-width: 0.75; fill: #D3D3D3;' /> <rect x='374.40' y='75.47' width='146.00' height='410.67' style='stroke-width: 0.75; fill: #D3D3D3;' /> <rect x='520.40' y='349.24' width='146.00' height='136.89' style='stroke-width: 0.75; fill: #D3D3D3;' /> -<line x1='2159.04' y1='502.56' x2='2159.04' y2='59.04' style='stroke-width: 0.75; stroke: #DF536B;' /> +<line x1='232.06' y1='502.56' x2='232.06' y2='59.04' style='stroke-width: 0.75; stroke: #DF536B;' /> <line x1='101.38' y1='95.62' x2='122.98' y2='95.62' style='stroke-width: 0.75; stroke: #DF536B;' /> <text x='133.78' y='99.75' style='font-size: 12.00px; font-family: sans;' textLength='51.35px' lengthAdjust='spacingAndGlyphs'>original fit</text> </g> diff --git a/tests/testthat/_snaps/multistart/parhist-for-biphasic-saemix-fit.svg b/tests/testthat/_snaps/multistart/parplot-for-biphasic-saemix-fit.svg index 75519b07..75519b07 100644 --- a/tests/testthat/_snaps/multistart/parhist-for-biphasic-saemix-fit.svg +++ b/tests/testthat/_snaps/multistart/parplot-for-biphasic-saemix-fit.svg diff --git a/tests/testthat/_snaps/multistart/parhist-for-sfo-fit.svg b/tests/testthat/_snaps/multistart/parplot-for-sfo-fit.svg index f3373901..f3373901 100644 --- a/tests/testthat/_snaps/multistart/parhist-for-sfo-fit.svg +++ b/tests/testthat/_snaps/multistart/parplot-for-sfo-fit.svg diff --git a/tests/testthat/test_multistart.R b/tests/testthat/test_multistart.R index e4885191..56eb140c 100644 --- a/tests/testthat/test_multistart.R +++ b/tests/testthat/test_multistart.R @@ -25,14 +25,14 @@ test_that("multistart works for saem.mmkin models", { skip_on_travis() # Plots are platform dependent llhist_sfo <- function() llhist(saem_sfo_s_multi) - parhist_sfo <- function() parhist(saem_sfo_s_multi, ylim = c(0.5, 2)) + parplot_sfo <- function() parplot(saem_sfo_s_multi, ylim = c(0.5, 2)) vdiffr::expect_doppelganger("llhist for sfo fit", llhist_sfo) - vdiffr::expect_doppelganger("parhist for sfo fit", parhist_sfo) + vdiffr::expect_doppelganger("parplot for sfo fit", parplot_sfo) llhist_biphasic <- function() llhist(saem_biphasic_m_multi) - parhist_biphasic <- function() parhist(saem_biphasic_m_multi, + parplot_biphasic <- function() parplot(saem_biphasic_m_multi, ylim = c(0.5, 2)) vdiffr::expect_doppelganger("llhist for biphasic saemix fit", llhist_biphasic) - vdiffr::expect_doppelganger("parhist for biphasic saemix fit", parhist_biphasic) + vdiffr::expect_doppelganger("parplot for biphasic saemix fit", parplot_biphasic) }) diff --git a/vignettes/web_only/multistart.rmd b/vignettes/web_only/multistart.rmd index 94f4ef98..27a8a96a 100644 --- a/vignettes/web_only/multistart.rmd +++ b/vignettes/web_only/multistart.rmd @@ -40,7 +40,7 @@ with different starting values. ```{r} f_saem_full_multi <- multistart(f_saem_full, n = 16, cores = 16) -parhist(f_saem_full_multi) +parplot(f_saem_full_multi) ``` This confirms that the variance of k2 is the most problematic parameter, so we @@ -51,7 +51,7 @@ for k2. f_saem_reduced <- update(f_saem_full, no_random_effect = "log_k2") illparms(f_saem_reduced) f_saem_reduced_multi <- multistart(f_saem_reduced, n = 16, cores = 16) -parhist(f_saem_reduced_multi, lpos = "topright") +parplot(f_saem_reduced_multi, lpos = "topright") ``` The results confirm that all remaining parameters can be determined with sufficient @@ -67,7 +67,7 @@ The parameter histograms can be further improved by excluding the result with the low likelihood. ```{r} -parhist(f_saem_reduced_multi, lpos = "topright", llmin = -326, ylim = c(0.5, 2)) +parplot(f_saem_reduced_multi, lpos = "topright", llmin = -326, ylim = c(0.5, 2)) ``` We can use the `anova` method to compare the models, including a likelihood ratio |
