aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--NAMESPACE4
-rw-r--r--R/mkinpredict.R10
-rw-r--r--R/saemix.R85
-rw-r--r--_pkgdown.yml1
-rw-r--r--docs/dev/pkgdown.yml2
-rw-r--r--docs/dev/sitemap.xml2
-rw-r--r--man/mkinpredict.Rd22
-rw-r--r--man/saem.Rd (renamed from man/saemix.Rd)77
-rw-r--r--test.log21
-rw-r--r--tests/testthat/FOCUS_2006_D.csf2
10 files changed, 103 insertions, 123 deletions
diff --git a/NAMESPACE b/NAMESPACE
index 33f3dec0..3360f861 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -28,7 +28,7 @@ S3method(print,nlme.mmkin)
S3method(print,summary.mkinfit)
S3method(print,summary.nlme.mmkin)
S3method(residuals,mkinfit)
-S3method(saemix,mmkin)
+S3method(saem,mmkin)
S3method(summary,mkinfit)
S3method(summary,nlme.mmkin)
S3method(update,mkinfit)
@@ -77,7 +77,7 @@ export(parms)
export(plot_err)
export(plot_res)
export(plot_sep)
-export(saemix)
+export(saem)
export(saemix_data)
export(saemix_model)
export(sigma_twocomp)
diff --git a/R/mkinpredict.R b/R/mkinpredict.R
index a6e7ca1c..7222e247 100644
--- a/R/mkinpredict.R
+++ b/R/mkinpredict.R
@@ -103,12 +103,7 @@
#' }
#'
#' @export
-mkinpredict <- function(x, odeparms, odeini,
- outtimes = seq(0, 120, by = 0.1),
- solution_type = "deSolve",
- use_compiled = "auto",
- method.ode = "lsoda", atol = 1e-8, rtol = 1e-10,
- map_output = TRUE, ...)
+mkinpredict <- function(x, odeparms, odeini, outtimes, ...)
{
UseMethod("mkinpredict", x)
}
@@ -122,8 +117,9 @@ mkinpredict.mkinmod <- function(x,
solution_type = "deSolve",
use_compiled = "auto",
method.ode = "lsoda", atol = 1e-8, rtol = 1e-10,
+ map_output = TRUE,
na_stop = TRUE,
- map_output = TRUE, ...)
+ ...)
{
# Names of state variables and observed variables
diff --git a/R/saemix.R b/R/saemix.R
index 280490a0..1db8b011 100644
--- a/R/saemix.R
+++ b/R/saemix.R
@@ -1,20 +1,16 @@
-#' Create saemix models
+#' Fit nonlinear mixed models with SAEM
#'
-#' The saemix function defined in this package is an S3 generic function
-#' using [saemix::saemix()] as its method for [saemix::SaemixModel] objects.
+#' This function uses [saemix::saemix()] as a backend for fitting nonlinear mixed
+#' effects models created from [mmkin] row objects using the stochastic approximation
+#' to the expectation maximisation algorithm (SAEM).
#'
-#' The method for mmkin row objects sets up a nonlinear mixed effects model for
-#' use with the saemix package. An mmkin row object is essentially a list of
-#' mkinfit objects that have been obtained by fitting the same model to
-#' a list of datasets.
+#' An mmkin row object is essentially a list of mkinfit objects that have been
+#' obtained by fitting the same model to a list of datasets using [mkinfit].
#'
-#' Starting values for the fixed effects (population mean parameters, argument psi0 of
-#' [saemix::saemixModel()] are the mean values of the parameters found using
-#' [mmkin].
+#' Starting values for the fixed effects (population mean parameters, argument
+#' psi0 of [saemix::saemixModel()] are the mean values of the parameters found
+#' using [mmkin].
#'
-#' @param model For the default method, this is an [saemix::saemixModel] object.
-#' If this is an [mmkin] row object, the [saemix::saemixModel] is created
-#' internally from the [mmkin] object.
#' @param object An [mmkin] row object containing several fits of the same
#' [mkinmod] model to different datasets
#' @param verbose Should we print information about created objects?
@@ -22,75 +18,84 @@
#' [parallel::mclapply()]. Using more than 1 core is experimental and may
#' lead to uncontrolled forking, apparently depending on the BLAS version
#' used.
+#' @param suppressPlot Should we suppress any plotting that is done
+#' by the saemix function?
+#' @param control Passed to [saemix::saemix]
#' @param \dots Further parameters passed to [saemix::saemixData]
#' and [saemix::saemixModel].
#' @return An [saemix::SaemixObject].
#' @examples
#' \dontrun{
-#' # We can load saemix, but should exclude the saemix function
-#' # as it would mask our generic version of it
-#' library(saemix, exclude = "saemix")
#' ds <- lapply(experimental_data_for_UBA_2019[6:10],
#' function(x) subset(x$data[c("name", "time", "value")]))
#' names(ds) <- paste("Dataset", 6:10)
#' f_mmkin_parent_p0_fixed <- mmkin("FOMC", ds, cores = 1,
#' state.ini = c(parent = 100), fixed_initials = "parent", quiet = TRUE)
-#' f_saemix_p0_fixed <- saemix(f_mmkin_parent_p0_fixed)
+#' f_saem_p0_fixed <- saem(f_mmkin_parent_p0_fixed)
#'
#' f_mmkin_parent <- mmkin(c("SFO", "FOMC", "DFOP"), ds, quiet = TRUE)
-#' f_saemix_sfo <- saemix(f_mmkin_parent["SFO", ])
-#' f_saemix_fomc <- saemix(f_mmkin_parent["FOMC", ])
-#' f_saemix_dfop <- saemix(f_mmkin_parent["DFOP", ])
+#' f_saem_sfo <- saem(f_mmkin_parent["SFO", ])
+#' f_saem_fomc <- saem(f_mmkin_parent["FOMC", ])
+#' f_saem_dfop <- saem(f_mmkin_parent["DFOP", ])
#'
-#' # As this returns an SaemixObject, we can use functions from saemix
-#' compare.saemix(list(f_saemix_sfo, f_saemix_fomc, f_saemix_dfop))
+#' # The returned saem.mmkin object contains an SaemixObject, we can use
+#' # functions from saemix
+#' library(saemix)
+#' compare.saemix(list(f_saem_sfo$so, f_saem_fomc$so, f_saem_dfop$so))
#'
#' f_mmkin_parent_tc <- update(f_mmkin_parent, error_model = "tc")
-#' f_saemix_fomc_tc <- saemix(f_mmkin_parent_tc["FOMC", ])
-#' compare.saemix(list(f_saemix_fomc, f_saemix_fomc_tc))
+#' f_saem_fomc_tc <- saem(f_mmkin_parent_tc["FOMC", ])
+#' compare.saemix(list(f_saem_fomc$so, f_saem_fomc_tc$so))
#'
#' dfop_sfo <- mkinmod(parent = mkinsub("DFOP", "A1"),
#' A1 = mkinsub("SFO"))
#' f_mmkin <- mmkin(list("DFOP-SFO" = dfop_sfo), ds, quiet = TRUE, solution_type = "analytical")
#' # This takes about 4 minutes on my system
-#' f_saemix <- saemix(f_mmkin)
+#' f_saem <- saem(f_mmkin)
#'
-#' # Using a single core, it takes about 6 minutes, using 10 cores it is slower
-#' # instead of faster
#' f_mmkin_des <- mmkin(list("DFOP-SFO" = dfop_sfo), ds, quiet = TRUE, solution_type = "deSolve")
-#' f_saemix_des <- saemix(f_mmkin_des, cores = 1)
-#' compare.saemix(list(f_saemix, f_saemix_des))
-#'
+#' # Using a single core, the following takes about 6 minutes, using 10 cores
+#' # it is slower instead of faster
+#' f_saem_des <- saem(f_mmkin_des, cores = 1)
+#' compare.saemix(list(f_saemix$so, f_saemix_des$so))
#' }
#' @export
-saemix <- function(model, data, control, ...) UseMethod("saemix")
+saem <- function(object, control, ...) UseMethod("saem")
-#' @rdname saemix
+#' @rdname saem
#' @export
-saemix.mmkin <- function(model, data,
+saem.mmkin <- function(object,
control = list(displayProgress = FALSE, print = FALSE,
save = FALSE, save.graphs = FALSE),
cores = 1,
verbose = FALSE, suppressPlot = TRUE, ...)
{
- m_saemix <- saemix_model(model, cores = cores, verbose = verbose)
- d_saemix <- saemix_data(model, verbose = verbose)
+ m_saemix <- saemix_model(object, cores = cores, verbose = verbose)
+ d_saemix <- saemix_data(object, verbose = verbose)
if (suppressPlot) {
# We suppress the log-likelihood curve that saemix currently
# produces at the end of the fit by plotting to a file
# that we discard afterwards
tmp <- tempfile()
- png(tmp)
+ grDevices::png(tmp)
}
- result <- saemix::saemix(m_saemix, d_saemix, control)
+ f_saemix <- saemix::saemix(m_saemix, d_saemix, control)
if (suppressPlot) {
- dev.off()
+ grDevices::dev.off()
unlink(tmp)
}
+ result <- list(
+ mkinmod = object[[1]]$mkinmod,
+ mmkin = object,
+ solution_type = object[[1]]$solution_type,
+ transform_rates = object[[1]]$transform_rates,
+ transform_fractions = object[[1]]$transform_fractions,
+ so = f_saemix)
+ class(result) <- "saem.mmkin"
return(result)
}
-#' @rdname saemix
+#' @rdname saem
#' @return An [saemix::SaemixModel] object.
#' @export
saemix_model <- function(object, cores = 1, verbose = FALSE, ...) {
@@ -254,7 +259,7 @@ saemix_model <- function(object, cores = 1, verbose = FALSE, ...) {
return(res)
}
-#' @rdname saemix
+#' @rdname saem
#' @return An [saemix::SaemixData] object.
#' @export
saemix_data <- function(object, verbose = FALSE, ...) {
diff --git a/_pkgdown.yml b/_pkgdown.yml
index ac100ac1..3668750f 100644
--- a/_pkgdown.yml
+++ b/_pkgdown.yml
@@ -44,6 +44,7 @@ reference:
- nlme.mmkin
- plot.nlme.mmkin
- summary.nlme.mmkin
+ - saem.mmkin
- nlme_function
- get_deg_func
- saemix_model
diff --git a/docs/dev/pkgdown.yml b/docs/dev/pkgdown.yml
index 207aaf78..e5410fd5 100644
--- a/docs/dev/pkgdown.yml
+++ b/docs/dev/pkgdown.yml
@@ -10,7 +10,7 @@ articles:
web_only/NAFTA_examples: NAFTA_examples.html
web_only/benchmarks: benchmarks.html
web_only/compiled_models: compiled_models.html
-last_built: 2020-11-05T23:08Z
+last_built: 2020-11-07T11:42Z
urls:
reference: https://pkgdown.jrwb.de/mkin/reference
article: https://pkgdown.jrwb.de/mkin/articles
diff --git a/docs/dev/sitemap.xml b/docs/dev/sitemap.xml
index 67e100de..c5f285f1 100644
--- a/docs/dev/sitemap.xml
+++ b/docs/dev/sitemap.xml
@@ -175,7 +175,7 @@
<loc>https://pkgdown.jrwb.de/mkin/reference/residuals.mkinfit.html</loc>
</url>
<url>
- <loc>https://pkgdown.jrwb.de/mkin/reference/saemix.html</loc>
+ <loc>https://pkgdown.jrwb.de/mkin/reference/saem.html</loc>
</url>
<url>
<loc>https://pkgdown.jrwb.de/mkin/reference/schaefer07_complex_case.html</loc>
diff --git a/man/mkinpredict.Rd b/man/mkinpredict.Rd
index f2799bf4..3a2939e1 100644
--- a/man/mkinpredict.Rd
+++ b/man/mkinpredict.Rd
@@ -6,19 +6,7 @@
\alias{mkinpredict.mkinfit}
\title{Produce predictions from a kinetic model using specific parameters}
\usage{
-mkinpredict(
- x,
- odeparms,
- odeini,
- outtimes = seq(0, 120, by = 0.1),
- solution_type = "deSolve",
- use_compiled = "auto",
- method.ode = "lsoda",
- atol = 1e-08,
- rtol = 1e-10,
- map_output = TRUE,
- ...
-)
+mkinpredict(x, odeparms, odeini, outtimes, ...)
\method{mkinpredict}{mkinmod}(
x,
@@ -30,8 +18,8 @@ mkinpredict(
method.ode = "lsoda",
atol = 1e-08,
rtol = 1e-10,
- na_stop = TRUE,
map_output = TRUE,
+ na_stop = TRUE,
...
)
@@ -65,6 +53,9 @@ observed variables, for example in the case of the SFORB model.}
\item{outtimes}{A numeric vector specifying the time points for which model
predictions should be generated.}
+\item{\dots}{Further arguments passed to the ode solver in case such a
+solver is used.}
+
\item{solution_type}{The method that should be used for producing the
predictions. This should generally be "analytical" if there is only one
observed variable, and usually "deSolve" in the case of several observed
@@ -89,9 +80,6 @@ the observed variables (default) or for all state variables (if set to
FALSE). Setting this to FALSE has no effect for analytical solutions,
as these always return mapped output.}
-\item{\dots}{Further arguments passed to the ode solver in case such a
-solver is used.}
-
\item{na_stop}{Should it be an error if deSolve::ode returns NaN values}
}
\value{
diff --git a/man/saemix.Rd b/man/saem.Rd
index a664b0cc..39b66448 100644
--- a/man/saemix.Rd
+++ b/man/saem.Rd
@@ -1,17 +1,16 @@
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/saemix.R
-\name{saemix}
-\alias{saemix}
-\alias{saemix.mmkin}
+\name{saem}
+\alias{saem}
+\alias{saem.mmkin}
\alias{saemix_model}
\alias{saemix_data}
-\title{Create saemix models}
+\title{Fit nonlinear mixed models with SAEM}
\usage{
-saemix(model, data, control, ...)
+saem(object, control, ...)
-\method{saemix}{mmkin}(
- model,
- data,
+\method{saem}{mmkin}(
+ object,
control = list(displayProgress = FALSE, print = FALSE, save = FALSE, save.graphs =
FALSE),
cores = 1,
@@ -25,9 +24,10 @@ saemix_model(object, cores = 1, verbose = FALSE, ...)
saemix_data(object, verbose = FALSE, ...)
}
\arguments{
-\item{model}{For the default method, this is an \link[saemix:saemixModel]{saemix::saemixModel} object.
-If this is an \link{mmkin} row object, the \link[saemix:saemixModel]{saemix::saemixModel} is created
-internally from the \link{mmkin} object.}
+\item{object}{An \link{mmkin} row object containing several fits of the same
+\link{mkinmod} model to different datasets}
+
+\item{control}{Passed to \link[saemix:saemix]{saemix::saemix}}
\item{\dots}{Further parameters passed to \link[saemix:saemixData]{saemix::saemixData}
and \link[saemix:saemixModel]{saemix::saemixModel}.}
@@ -39,8 +39,8 @@ used.}
\item{verbose}{Should we print information about created objects?}
-\item{object}{An \link{mmkin} row object containing several fits of the same
-\link{mkinmod} model to different datasets}
+\item{suppressPlot}{Should we suppress any plotting that is done
+by the saemix function?}
}
\value{
An \link[saemix:SaemixObject-class]{saemix::SaemixObject}.
@@ -50,54 +50,51 @@ An \link[saemix:SaemixModel-class]{saemix::SaemixModel} object.
An \link[saemix:SaemixData-class]{saemix::SaemixData} object.
}
\description{
-The saemix function defined in this package is an S3 generic function
-using \code{\link[saemix:saemix]{saemix::saemix()}} as its method for \link[saemix:SaemixModel-class]{saemix::SaemixModel} objects.
+This function uses \code{\link[saemix:saemix]{saemix::saemix()}} as a backend for fitting nonlinear mixed
+effects models created from \link{mmkin} row objects using the stochastic approximation
+to the expectation maximisation algorithm (SAEM).
}
\details{
-The method for mmkin row objects sets up a nonlinear mixed effects model for
-use with the saemix package. An mmkin row object is essentially a list of
-mkinfit objects that have been obtained by fitting the same model to
-a list of datasets.
-
-Starting values for the fixed effects (population mean parameters, argument psi0 of
-\code{\link[saemix:saemixModel]{saemix::saemixModel()}} are the mean values of the parameters found using
-\link{mmkin}.
+An mmkin row object is essentially a list of mkinfit objects that have been
+obtained by fitting the same model to a list of datasets using \link{mkinfit}.
+
+Starting values for the fixed effects (population mean parameters, argument
+psi0 of \code{\link[saemix:saemixModel]{saemix::saemixModel()}} are the mean values of the parameters found
+using \link{mmkin}.
}
\examples{
\dontrun{
-# We can load saemix, but should exclude the saemix function
-# as it would mask our generic version of it
-library(saemix, exclude = "saemix")
ds <- lapply(experimental_data_for_UBA_2019[6:10],
function(x) subset(x$data[c("name", "time", "value")]))
names(ds) <- paste("Dataset", 6:10)
f_mmkin_parent_p0_fixed <- mmkin("FOMC", ds, cores = 1,
state.ini = c(parent = 100), fixed_initials = "parent", quiet = TRUE)
-f_saemix_p0_fixed <- saemix(f_mmkin_parent_p0_fixed)
+f_saem_p0_fixed <- saem(f_mmkin_parent_p0_fixed)
f_mmkin_parent <- mmkin(c("SFO", "FOMC", "DFOP"), ds, quiet = TRUE)
-f_saemix_sfo <- saemix(f_mmkin_parent["SFO", ])
-f_saemix_fomc <- saemix(f_mmkin_parent["FOMC", ])
-f_saemix_dfop <- saemix(f_mmkin_parent["DFOP", ])
+f_saem_sfo <- saem(f_mmkin_parent["SFO", ])
+f_saem_fomc <- saem(f_mmkin_parent["FOMC", ])
+f_saem_dfop <- saem(f_mmkin_parent["DFOP", ])
-# As this returns an SaemixObject, we can use functions from saemix
-compare.saemix(list(f_saemix_sfo, f_saemix_fomc, f_saemix_dfop))
+# The returned saem.mmkin object contains an SaemixObject, we can use
+# functions from saemix
+library(saemix)
+compare.saemix(list(f_saem_sfo$so, f_saem_fomc$so, f_saem_dfop$so))
f_mmkin_parent_tc <- update(f_mmkin_parent, error_model = "tc")
-f_saemix_fomc_tc <- saemix(f_mmkin_parent_tc["FOMC", ])
-compare.saemix(list(f_saemix_fomc, f_saemix_fomc_tc))
+f_saem_fomc_tc <- saem(f_mmkin_parent_tc["FOMC", ])
+compare.saemix(list(f_saem_fomc$so, f_saem_fomc_tc$so))
dfop_sfo <- mkinmod(parent = mkinsub("DFOP", "A1"),
A1 = mkinsub("SFO"))
f_mmkin <- mmkin(list("DFOP-SFO" = dfop_sfo), ds, quiet = TRUE, solution_type = "analytical")
# This takes about 4 minutes on my system
-f_saemix <- saemix(f_mmkin)
+f_saem <- saem(f_mmkin)
-# Using a single core, it takes about 6 minutes, using 10 cores it is slower
-# instead of faster
f_mmkin_des <- mmkin(list("DFOP-SFO" = dfop_sfo), ds, quiet = TRUE, solution_type = "deSolve")
-f_saemix_des <- saemix(f_mmkin_des, cores = 1)
-compare.saemix(list(f_saemix, f_saemix_des))
-
+# Using a single core, the following takes about 6 minutes, using 10 cores
+# it is slower instead of faster
+f_saem_des <- saem(f_mmkin_des, cores = 1)
+compare.saemix(list(f_saemix$so, f_saemix_des$so))
}
}
diff --git a/test.log b/test.log
index 318cc304..22501b25 100644
--- a/test.log
+++ b/test.log
@@ -5,11 +5,11 @@ Testing mkin
✔ | 2 | Export dataset for reading into CAKE
✔ | 14 | Results for FOCUS D established in expertise for UBA (Ranke 2014) [0.9 s]
✔ | 4 | Calculation of FOCUS chi2 error levels [0.4 s]
-✔ | 7 | Fitting the SFORB model [3.4 s]
-✔ | 5 | Analytical solutions for coupled models [3.1 s]
+✔ | 7 | Fitting the SFORB model [3.3 s]
+✔ | 5 | Analytical solutions for coupled models [3.0 s]
✔ | 5 | Calculation of Akaike weights
✔ | 10 | Confidence intervals and p-values [1.0 s]
-✔ | 14 | Error model fitting [4.3 s]
+✔ | 14 | Error model fitting [4.1 s]
✔ | 4 | Test fitting the decline of metabolites from their maximum [0.3 s]
✔ | 1 | Fitting the logistic model [0.2 s]
✔ | 1 | Test dataset class mkinds used in gmkin
@@ -33,23 +33,16 @@ Reason: getRversion() < "4.1.0" is TRUE
────────────────────────────────────────────────────────────────────────────────
✔ | 4 | Residuals extracted from mkinfit models
✔ | 2 | Complex test case from Schaefer et al. (2007) Piacenza paper [1.5 s]
-✖ | 3 1 | Summary [0.1 s]
-────────────────────────────────────────────────────────────────────────────────
-test_summary.R:17: failure: Summaries are reproducible
-Results have changed from known value recorded in 'summary_DFOP_FOCUS_C.txt'.
-1/82 mismatches
-x[66]: "parent 2.661 4 5"
-y[66]: "parent 2.495 3 6"
-────────────────────────────────────────────────────────────────────────────────
+✔ | 4 | Summary [0.1 s]
✔ | 1 | Summaries of old mkinfit objects
✔ | 4 | Results for synthetic data established in expertise for UBA (Ranke 2014) [2.2 s]
✔ | 9 | Hypothesis tests [7.0 s]
✔ | 4 | Calculation of maximum time weighted average concentrations (TWAs) [2.5 s]
══ Results ═════════════════════════════════════════════════════════════════════
-Duration: 38.1 s
+Duration: 37.7 s
-OK: 145
-Failed: 1
+OK: 146
+Failed: 0
Warnings: 0
Skipped: 3
diff --git a/tests/testthat/FOCUS_2006_D.csf b/tests/testthat/FOCUS_2006_D.csf
index 8da4cd91..dc369ee1 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-11-05
+Date: 2020-11-07
Optimiser: IRLS
[Data]

Contact - Imprint