aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorJohannes Ranke <jranke@uni-bremen.de>2022-10-28 13:39:15 +0200
committerJohannes Ranke <jranke@uni-bremen.de>2022-10-28 13:39:15 +0200
commitf820bf5b91be0f589de16c3e3250f5f79672df75 (patch)
tree2b1406e1c9286634ca017db586e09e2299dec048
parentb1740ade9a1746ccdb325b95915ef88872489f03 (diff)
Rename parhist to parplot and make it generic
That parhist name was not the brightest idea, as it does not show histograms.
-rw-r--r--NAMESPACE3
-rw-r--r--NEWS.md4
-rw-r--r--R/llhist.R3
-rw-r--r--R/multistart.R6
-rw-r--r--R/parplot.R (renamed from R/parhist.R)12
-rw-r--r--R/saem.R2
-rw-r--r--log/test.log42
-rw-r--r--man/logLik.saem.mmkin.Rd2
-rw-r--r--man/multistart.Rd6
-rw-r--r--man/parplot.Rd (renamed from man/parhist.Rd)17
-rw-r--r--tests/testthat/_snaps/multistart/llhist-for-biphasic-saemix-fit.svg1
-rw-r--r--tests/testthat/_snaps/multistart/llhist-for-sfo-fit.svg2
-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.R8
-rw-r--r--vignettes/web_only/multistart.rmd6
16 files changed, 63 insertions, 51 deletions
diff --git a/NAMESPACE b/NAMESPACE
index a5691eee..4755f631 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -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)
diff --git a/NEWS.md b/NEWS.md
index d288c32f..c7927559 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -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.
diff --git a/R/llhist.R b/R/llhist.R
index 22e3aa08..e158495d 100644
--- a/R/llhist.R
+++ b/R/llhist.R
@@ -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)
diff --git a/R/saem.R b/R/saem.R
index cf67b8e1..090ed3bf 100644
--- a/R/saem.R
+++ b/R/saem.R
@@ -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

Contact - Imprint