From 48c463680b51fa767b4cd7bd62865f192d0354ac Mon Sep 17 00:00:00 2001
From: Johannes Ranke
Date: Sat, 6 Feb 2021 18:30:32 +0100
Subject: Reintroduce interface to saemix
Also after the upgrade from buster to bullseye of my local system, some
test results for saemix have changed.
---
.travis.yml | 2 -
DESCRIPTION | 5 +-
NAMESPACE | 8 +
NEWS.md | 10 +
R/endpoints.R | 8 +-
R/plot.mixed.mmkin.R | 23 +-
R/saem.R | 512 +++
R/summary.saem.mmkin.R | 268 ++
_pkgdown.yml | 7 +-
build.log | 2 +-
check.log | 15 +-
docs/dev/404.html | 2 +-
docs/dev/articles/index.html | 2 +-
docs/dev/authors.html | 8 +-
docs/dev/index.html | 4 +-
docs/dev/news/index.html | 211 +-
docs/dev/pkgdown.yml | 4 +-
docs/dev/reference/Rplot001.png | Bin 19395 -> 1011 bytes
docs/dev/reference/Rplot003.png | Bin 28735 -> 28745 bytes
docs/dev/reference/Rplot005.png | Bin 57095 -> 59143 bytes
docs/dev/reference/Rplot006.png | Bin 24624 -> 22103 bytes
docs/dev/reference/confint.mkinfit.html | 127 +-
docs/dev/reference/endpoints.html | 4 +-
docs/dev/reference/index.html | 10 +-
docs/dev/reference/logLik.mkinfit.html | 12 +-
docs/dev/reference/mkinresplot.html | 6 +-
docs/dev/reference/mmkin-1.png | Bin 114048 -> 110459 bytes
docs/dev/reference/mmkin-2.png | Bin 110392 -> 107057 bytes
docs/dev/reference/mmkin-3.png | Bin 97556 -> 96062 bytes
docs/dev/reference/mmkin-4.png | Bin 70005 -> 67191 bytes
docs/dev/reference/mmkin-5.png | Bin 66093 -> 64880 bytes
docs/dev/reference/mmkin.html | 31 +-
docs/dev/reference/nlme.mmkin-1.png | Bin 119655 -> 124677 bytes
docs/dev/reference/nlme.mmkin-2.png | Bin 159253 -> 169523 bytes
docs/dev/reference/nlme.mmkin-3.png | Bin 161800 -> 172692 bytes
docs/dev/reference/nlme.mmkin.html | 76 +-
docs/dev/reference/plot.mixed.mmkin-1.png | Bin 86076 -> 84734 bytes
docs/dev/reference/plot.mixed.mmkin-2.png | Bin 164488 -> 173916 bytes
docs/dev/reference/plot.mixed.mmkin-3.png | Bin 163536 -> 172396 bytes
docs/dev/reference/plot.mixed.mmkin-4.png | Bin 166687 -> 175502 bytes
docs/dev/reference/plot.mixed.mmkin.html | 8 +-
docs/dev/reference/saem-3.png | Bin 82102 -> 82107 bytes
docs/dev/reference/saem-5.png | Bin 164443 -> 173288 bytes
docs/dev/reference/saem.html | 56 +-
docs/dev/reference/summary.saem.mmkin.html | 426 +-
docs/dev/reference/transform_odeparms.html | 80 +-
docs/dev/sitemap.xml | 3 -
man/endpoints.Rd | 8 +-
man/plot.mixed.mmkin.Rd | 11 +-
man/saem.Rd | 155 +
man/summary.saem.mmkin.Rd | 100 +
test.log | 27 +-
...t-for-saem-object-with-mkin-transformations.svg | 4782 ++++++++++----------
...for-saem-object-with-saemix-transformations.svg | 1234 ++---
tests/testthat/setup_script.R | 18 +
tests/testthat/summary_saem_biphasic_s.txt | 77 +
tests/testthat/test_mixed.R | 135 +-
tests/testthat/test_plot.R | 11 +
58 files changed, 4938 insertions(+), 3550 deletions(-)
create mode 100644 R/saem.R
create mode 100644 R/summary.saem.mmkin.R
create mode 100644 man/saem.Rd
create mode 100644 man/summary.saem.mmkin.Rd
create mode 100644 tests/testthat/summary_saem_biphasic_s.txt
diff --git a/.travis.yml b/.travis.yml
index 282f1302..6c03b451 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -11,8 +11,6 @@ addons:
cache: packages
repos:
CRAN: https://cloud.r-project.org
-r_github_packages:
- - saemixdevelopment/saemixextension@master
script:
- R CMD build .
- R CMD check --no-tests mkin_*.tar.gz
diff --git a/DESCRIPTION b/DESCRIPTION
index dd57aefa..15871070 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
Package: mkin
Type: Package
Title: Kinetic Evaluation of Chemical Degradation Data
-Version: 1.0.1
+Version: 1.0.1.9000
Date: 2021-02-06
Authors@R: c(
person("Johannes", "Ranke", role = c("aut", "cre", "cph"),
@@ -19,9 +19,10 @@ Description: Calculation routines based on the FOCUS Kinetics Report (2006,
particular purpose.
Depends: R (>= 2.15.1), parallel
Imports: stats, graphics, methods, deSolve, R6, inline (>= 0.3.17), numDeriv,
- lmtest, pkgbuild, nlme (>= 3.1-151), purrr
+ lmtest, pkgbuild, nlme (>= 3.1-151), purrr, saemix (>= 3.1.9000)
Suggests: knitr, rbenchmark, tikzDevice, testthat, rmarkdown, covr, vdiffr,
benchmarkme, tibble, stats4
+Remotes: github::saemixdevelopment/saemixextension
License: GPL
LazyLoad: yes
LazyData: yes
diff --git a/NAMESPACE b/NAMESPACE
index 9776d2f3..f2497283 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -30,11 +30,15 @@ S3method(print,mkinmod)
S3method(print,mmkin)
S3method(print,nafta)
S3method(print,nlme.mmkin)
+S3method(print,saem.mmkin)
S3method(print,summary.mkinfit)
S3method(print,summary.nlme.mmkin)
+S3method(print,summary.saem.mmkin)
S3method(residuals,mkinfit)
+S3method(saem,mmkin)
S3method(summary,mkinfit)
S3method(summary,nlme.mmkin)
+S3method(summary,saem.mmkin)
S3method(update,mkinfit)
S3method(update,mmkin)
S3method(update,nlme.mmkin)
@@ -86,6 +90,9 @@ export(parms)
export(plot_err)
export(plot_res)
export(plot_sep)
+export(saem)
+export(saemix_data)
+export(saemix_model)
export(sigma_twocomp)
export(transform_odeparms)
import(deSolve)
@@ -126,5 +133,6 @@ importFrom(stats,residuals)
importFrom(stats,rnorm)
importFrom(stats,shapiro.test)
importFrom(stats,update)
+importFrom(stats,vcov)
importFrom(utils,getFromNamespace)
importFrom(utils,write.table)
diff --git a/NEWS.md b/NEWS.md
index 124111d6..08745e9f 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,3 +1,13 @@
+# mkin 1.0.1.9000
+
+- Switch to a versioning scheme where the fourth version component indicates development versions
+
+- Reintroduce the interface to the current development version of saemix, in particular:
+
+- 'saemix_model' and 'saemix_data': Helper functions to set up nonlinear mixed-effects models for mmkin row objects
+
+- 'saem': generic function to fit saemix models using 'saemix_model' and 'saemix_data', with a generator 'saem.mmkin', summary and plot methods
+
# mkin 1.0.1
- 'confint.mmkin', 'nlme.mmkin', 'transform_odeparms': Fix example code in dontrun sections that failed with current defaults
diff --git a/R/endpoints.R b/R/endpoints.R
index b5872e68..f1f47581 100644
--- a/R/endpoints.R
+++ b/R/endpoints.R
@@ -10,8 +10,8 @@
#' Additional DT50 values are calculated from the FOMC DT90 and k1 and k2 from
#' HS and DFOP, as well as from Eigenvalues b1 and b2 of any SFORB models
#'
-#' @param fit An object of class [mkinfit] or [nlme.mmkin]
-#' or another object that has list components
+#' @param fit An object of class [mkinfit], [nlme.mmkin] or
+#' [saem.mmkin]. Or another object that has list components
#' mkinmod containing an [mkinmod] degradation model, and two numeric vectors,
#' bparms.optim and bparms.fixed, that contain parameter values
#' for that model.
@@ -20,8 +20,8 @@
#' and, if applicable, a vector of formation fractions named ff
#' and, if the SFORB model was in use, a vector of eigenvalues
#' of these SFORB models, equivalent to DFOP rate constants
-#' @note The function is used internally by [summary.mkinfit]
-#' and [summary.nlme.mmkin]
+#' @note The function is used internally by [summary.mkinfit],
+#' [summary.nlme.mmkin] and [summary.saem.mmkin].
#' @author Johannes Ranke
#' @examples
#'
diff --git a/R/plot.mixed.mmkin.R b/R/plot.mixed.mmkin.R
index 5a0b7412..1674d855 100644
--- a/R/plot.mixed.mmkin.R
+++ b/R/plot.mixed.mmkin.R
@@ -2,7 +2,7 @@ utils::globalVariables("ds")
#' Plot predictions from a fitted nonlinear mixed model obtained via an mmkin row object
#'
-#' @param x An object of class [mixed.mmkin], [nlme.mmkin]
+#' @param x An object of class [mixed.mmkin], [saem.mmkin] or [nlme.mmkin]
#' @param i A numeric index to select datasets for which to plot the individual predictions,
#' in case plots get too large
#' @inheritParams plot.mkinfit
@@ -39,6 +39,15 @@ utils::globalVariables("ds")
#' f_nlme <- nlme(f, control = list(pnlsMaxIter = 120, tolerance = 1e-3))
#' plot(f_nlme)
#'
+#' f_saem <- saem(f, transformations = "saemix")
+#' plot(f_saem)
+#'
+#' # We can overlay the two variants if we generate predictions
+#' pred_nlme <- mkinpredict(dfop_sfo,
+#' f_nlme$bparms.optim[-1],
+#' c(parent = f_nlme$bparms.optim[[1]], A1 = 0),
+#' seq(0, 180, by = 0.2))
+#' plot(f_saem, pred_over = list(nlme = pred_nlme))
#' }
#' @export
plot.mixed.mmkin <- function(x,
@@ -82,6 +91,18 @@ plot.mixed.mmkin <- function(x,
type = ifelse(standardized, "pearson", "response"))
}
+ if (inherits(x, "saem.mmkin")) {
+ if (x$transformations == "saemix") backtransform = FALSE
+ degparms_i <- saemix::psi(x$so)
+ rownames(degparms_i) <- ds_names
+ degparms_i_names <- setdiff(x$so@results@name.fixed, names(fit_1$errparms))
+ colnames(degparms_i) <- degparms_i_names
+ residual_type = ifelse(standardized, "standardized", "residual")
+ residuals <- x$data[[residual_type]]
+ degparms_pop <- x$so@results@fixed.effects
+ names(degparms_pop) <- degparms_i_names
+ }
+
degparms_fixed <- fit_1$fixed$value
names(degparms_fixed) <- rownames(fit_1$fixed)
degparms_all <- cbind(as.matrix(degparms_i),
diff --git a/R/saem.R b/R/saem.R
new file mode 100644
index 00000000..fd2a77b4
--- /dev/null
+++ b/R/saem.R
@@ -0,0 +1,512 @@
+utils::globalVariables(c("predicted", "std"))
+
+#' Fit nonlinear mixed models with SAEM
+#'
+#' This function uses [saemix::saemix()] as a backend for fitting nonlinear mixed
+#' effects models created from [mmkin] row objects using the Stochastic Approximation
+#' Expectation Maximisation algorithm (SAEM).
+#'
+#' 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].
+#'
+#' @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 of
+#' type [saemix::SaemixModel] and [saemix::SaemixData]?
+#' @param transformations Per default, all parameter transformations are done
+#' in mkin. If this argument is set to 'saemix', parameter transformations
+#' are done in 'saemix' for the supported cases. Currently this is only
+#' supported in cases where the initial concentration of the parent is not fixed,
+#' SFO or DFOP is used for the parent and there is either no metabolite or one.
+#' @param degparms_start Parameter values given as a named numeric vector will
+#' be used to override the starting values obtained from the 'mmkin' object.
+#' @param solution_type Possibility to specify the solution type in case the
+#' automatic choice is not desired
+#' @param quiet Should we suppress the messages saemix prints at the beginning
+#' and the end of the optimisation process?
+#' @param control Passed to [saemix::saemix]
+#' @param \dots Further parameters passed to [saemix::saemixModel].
+#' @return An S3 object of class 'saem.mmkin', containing the fitted
+#' [saemix::SaemixObject] as a list component named 'so'. The
+#' object also inherits from 'mixed.mmkin'.
+#' @seealso [summary.saem.mmkin] [plot.mixed.mmkin]
+#' @examples
+#' \dontrun{
+#' 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,
+#' state.ini = c(parent = 100), fixed_initials = "parent", quiet = TRUE)
+#' f_saem_p0_fixed <- saem(f_mmkin_parent_p0_fixed)
+#'
+#' f_mmkin_parent <- mmkin(c("SFO", "FOMC", "DFOP"), ds, quiet = TRUE)
+#' f_saem_sfo <- saem(f_mmkin_parent["SFO", ])
+#' f_saem_fomc <- saem(f_mmkin_parent["FOMC", ])
+#' f_saem_dfop <- saem(f_mmkin_parent["DFOP", ])
+#'
+#' # The returned saem.mmkin object contains an SaemixObject, therefore we can use
+#' # functions from saemix
+#' library(saemix)
+#' compare.saemix(list(f_saem_sfo$so, f_saem_fomc$so, f_saem_dfop$so))
+#' plot(f_saem_fomc$so, plot.type = "convergence")
+#' plot(f_saem_fomc$so, plot.type = "individual.fit")
+#' plot(f_saem_fomc$so, plot.type = "npde")
+#' plot(f_saem_fomc$so, plot.type = "vpc")
+#'
+#' f_mmkin_parent_tc <- update(f_mmkin_parent, error_model = "tc")
+#' f_saem_fomc_tc <- saem(f_mmkin_parent_tc["FOMC", ])
+#' compare.saemix(list(f_saem_fomc$so, f_saem_fomc_tc$so))
+#'
+#' sfo_sfo <- mkinmod(parent = mkinsub("SFO", "A1"),
+#' A1 = mkinsub("SFO"))
+#' fomc_sfo <- mkinmod(parent = mkinsub("FOMC", "A1"),
+#' A1 = mkinsub("SFO"))
+#' dfop_sfo <- mkinmod(parent = mkinsub("DFOP", "A1"),
+#' A1 = mkinsub("SFO"))
+#' # The following fit uses analytical solutions for SFO-SFO and DFOP-SFO,
+#' # and compiled ODEs for FOMC that are much slower
+#' f_mmkin <- mmkin(list(
+#' "SFO-SFO" = sfo_sfo, "FOMC-SFO" = fomc_sfo, "DFOP-SFO" = dfop_sfo),
+#' ds, quiet = TRUE)
+#' # saem fits of SFO-SFO and DFOP-SFO to these data take about five seconds
+#' # each on this system, as we use analytical solutions written for saemix.
+#' # When using the analytical solutions written for mkin this took around
+#' # four minutes
+#' f_saem_sfo_sfo <- saem(f_mmkin["SFO-SFO", ])
+#' f_saem_dfop_sfo <- saem(f_mmkin["DFOP-SFO", ])
+#' # We can use print, plot and summary methods to check the results
+#' print(f_saem_dfop_sfo)
+#' plot(f_saem_dfop_sfo)
+#' summary(f_saem_dfop_sfo, data = TRUE)
+#'
+#' # The following takes about 6 minutes
+#' #f_saem_dfop_sfo_deSolve <- saem(f_mmkin["DFOP-SFO", ], solution_type = "deSolve",
+#' # control = list(nbiter.saemix = c(200, 80), nbdisplay = 10))
+#'
+#' #saemix::compare.saemix(list(
+#' # f_saem_dfop_sfo$so,
+#' # f_saem_dfop_sfo_deSolve$so))
+#'
+#' # If the model supports it, we can also use eigenvalue based solutions, which
+#' # take a similar amount of time
+#' #f_saem_sfo_sfo_eigen <- saem(f_mmkin["SFO-SFO", ], solution_type = "eigen",
+#' # control = list(nbiter.saemix = c(200, 80), nbdisplay = 10))
+#' }
+#' @export
+saem <- function(object, ...) UseMethod("saem")
+
+#' @rdname saem
+#' @export
+saem.mmkin <- function(object,
+ transformations = c("mkin", "saemix"),
+ degparms_start = numeric(),
+ solution_type = "auto",
+ control = list(displayProgress = FALSE, print = FALSE,
+ save = FALSE, save.graphs = FALSE),
+ verbose = FALSE, quiet = FALSE, ...)
+{
+ transformations <- match.arg(transformations)
+ m_saemix <- saemix_model(object, verbose = verbose,
+ degparms_start = degparms_start, solution_type = solution_type,
+ transformations = transformations, ...)
+ d_saemix <- saemix_data(object, verbose = verbose)
+
+ fit_time <- system.time({
+ utils::capture.output(f_saemix <- saemix::saemix(m_saemix, d_saemix, control), split = !quiet)
+ })
+
+ transparms_optim <- f_saemix@results@fixed.effects
+ names(transparms_optim) <- f_saemix@results@name.fixed
+
+ if (transformations == "mkin") {
+ bparms_optim <- backtransform_odeparms(transparms_optim,
+ object[[1]]$mkinmod,
+ object[[1]]$transform_rates,
+ object[[1]]$transform_fractions)
+ } else {
+ bparms_optim <- transparms_optim
+ }
+
+ return_data <- nlme_data(object)
+
+ return_data$predicted <- f_saemix@model@model(
+ psi = saemix::psi(f_saemix),
+ id = as.numeric(return_data$ds),
+ xidep = return_data[c("time", "name")])
+
+ return_data <- transform(return_data,
+ residual = predicted - value,
+ std = sigma_twocomp(predicted,
+ f_saemix@results@respar[1], f_saemix@results@respar[2]))
+ return_data <- transform(return_data,
+ standardized = residual / std)
+
+ result <- list(
+ mkinmod = object[[1]]$mkinmod,
+ mmkin = object,
+ solution_type = object[[1]]$solution_type,
+ transformations = transformations,
+ transform_rates = object[[1]]$transform_rates,
+ transform_fractions = object[[1]]$transform_fractions,
+ so = f_saemix,
+ time = fit_time,
+ mean_dp_start = attr(m_saemix, "mean_dp_start"),
+ bparms.optim = bparms_optim,
+ bparms.fixed = object[[1]]$bparms.fixed,
+ data = return_data,
+ err_mod = object[[1]]$err_mod,
+ date.fit = date(),
+ saemixversion = as.character(utils::packageVersion("saemix")),
+ mkinversion = as.character(utils::packageVersion("mkin")),
+ Rversion = paste(R.version$major, R.version$minor, sep=".")
+ )
+
+ class(result) <- c("saem.mmkin", "mixed.mmkin")
+ return(result)
+}
+
+#' @export
+#' @rdname saem
+#' @param x An saem.mmkin object to print
+#' @param digits Number of digits to use for printing
+print.saem.mmkin <- function(x, digits = max(3, getOption("digits") - 3), ...) {
+ cat( "Kinetic nonlinear mixed-effects model fit by SAEM" )
+ cat("\nStructural model:\n")
+ diffs <- x$mmkin[[1]]$mkinmod$diffs
+ nice_diffs <- gsub("^(d.*) =", "\\1/dt =", diffs)
+ writeLines(strwrap(nice_diffs, exdent = 11))
+ cat("\nData:\n")
+ cat(nrow(x$data), "observations of",
+ length(unique(x$data$name)), "variable(s) grouped in",
+ length(unique(x$data$ds)), "datasets\n")
+
+ cat("\nLikelihood computed by importance sampling\n")
+ print(data.frame(
+ AIC = AIC(x$so, type = "is"),
+ BIC = BIC(x$so, type = "is"),
+ logLik = logLik(x$so, type = "is"),
+ row.names = " "), digits = digits)
+
+ cat("\nFitted parameters:\n")
+ conf.int <- x$so@results@conf.int[c("estimate", "lower", "upper")]
+ rownames(conf.int) <- x$so@results@conf.int[["name"]]
+ print(conf.int, digits = digits)
+
+ invisible(x)
+}
+
+#' @rdname saem
+#' @return An [saemix::SaemixModel] object.
+#' @export
+saemix_model <- function(object, solution_type = "auto", transformations = c("mkin", "saemix"),
+ degparms_start = numeric(), verbose = FALSE, ...)
+{
+ if (nrow(object) > 1) stop("Only row objects allowed")
+
+ mkin_model <- object[[1]]$mkinmod
+
+ degparms_optim <- mean_degparms(object)
+ if (transformations == "saemix") {
+ degparms_optim <- backtransform_odeparms(degparms_optim,
+ object[[1]]$mkinmod,
+ object[[1]]$transform_rates,
+ object[[1]]$transform_fractions)
+ }
+ degparms_fixed <- object[[1]]$bparms.fixed
+
+ # Transformations are done in the degradation function
+ transform.par = rep(0, length(degparms_optim))
+
+ odeini_optim_parm_names <- grep('_0$', names(degparms_optim), value = TRUE)
+ odeini_fixed_parm_names <- grep('_0$', names(degparms_fixed), value = TRUE)
+
+ odeparms_fixed_names <- setdiff(names(degparms_fixed), odeini_fixed_parm_names)
+ odeparms_fixed <- degparms_fixed[odeparms_fixed_names]
+
+ odeini_fixed <- degparms_fixed[odeini_fixed_parm_names]
+ names(odeini_fixed) <- gsub('_0$', '', odeini_fixed_parm_names)
+
+ model_function <- FALSE
+
+ # Model functions with analytical solutions
+ # Fixed parameters, use_of_ff = "min" and turning off sinks currently not supported here
+ # In general, we need to consider exactly how the parameters in mkinfit were specified,
+ # as the parameters are currently mapped by position in these solutions
+ sinks <- sapply(mkin_model$spec, function(x) x$sink)
+ if (length(odeparms_fixed) == 0 & mkin_model$use_of_ff == "max" & all(sinks)) {
+ # Parent only
+ if (length(mkin_model$spec) == 1) {
+ parent_type <- mkin_model$spec[[1]]$type
+ if (length(odeini_fixed) == 1) {
+ if (parent_type == "SFO") {
+ stop("saemix needs at least two parameters to work on.")
+ }
+ if (parent_type == "FOMC") {
+ model_function <- function(psi, id, xidep) {
+ odeini_fixed / (xidep[, "time"]/exp(psi[id, 2]) + 1)^exp(psi[id, 1])
+ }
+ }
+ if (parent_type == "DFOP") {
+ model_function <- function(psi, id, xidep) {
+ g <- plogis(psi[id, 3])
+ t <- xidep[, "time"]
+ odeini_fixed * (g * exp(- exp(psi[id, 1]) * t) +
+ (1 - g) * exp(- exp(psi[id, 2]) * t))
+ }
+ }
+ if (parent_type == "HS") {
+ model_function <- function(psi, id, xidep) {
+ tb <- exp(psi[id, 3])
+ t <- xidep[, "time"]
+ k1 = exp(psi[id, 1])
+ odeini_fixed * ifelse(t <= tb,
+ exp(- k1 * t),
+ exp(- k1 * tb) * exp(- exp(psi[id, 2]) * (t - tb)))
+ }
+ }
+ } else {
+ if (parent_type == "SFO") {
+ if (transformations == "mkin") {
+ model_function <- function(psi, id, xidep) {
+ psi[id, 1] * exp( - exp(psi[id, 2]) * xidep[, "time"])
+ }
+ } else {
+ model_function <- function(psi, id, xidep) {
+ psi[id, 1] * exp( - psi[id, 2] * xidep[, "time"])
+ }
+ transform.par = c(0, 1)
+ }
+ }
+ if (parent_type == "FOMC") {
+ model_function <- function(psi, id, xidep) {
+ psi[id, 1] / (xidep[, "time"]/exp(psi[id, 3]) + 1)^exp(psi[id, 2])
+ }
+ }
+ if (parent_type == "DFOP") {
+ if (transformations == "mkin") {
+ model_function <- function(psi, id, xidep) {
+ g <- plogis(psi[id, 4])
+ t <- xidep[, "time"]
+ psi[id, 1] * (g * exp(- exp(psi[id, 2]) * t) +
+ (1 - g) * exp(- exp(psi[id, 3]) * t))
+ }
+ } else {
+ model_function <- function(psi, id, xidep) {
+ g <- psi[id, 4]
+ t <- xidep[, "time"]
+ psi[id, 1] * (g * exp(- psi[id, 2] * t) +
+ (1 - g) * exp(- psi[id, 3] * t))
+ }
+ transform.par = c(0, 1, 1, 3)
+ }
+ }
+ if (parent_type == "HS") {
+ model_function <- function(psi, id, xidep) {
+ tb <- exp(psi[id, 4])
+ t <- xidep[, "time"]
+ k1 = exp(psi[id, 2])
+ psi[id, 1] * ifelse(t <= tb,
+ exp(- k1 * t),
+ exp(- k1 * tb) * exp(- exp(psi[id, 3]) * (t - tb)))
+ }
+ }
+ }
+ }
+
+ # Parent with one metabolite
+ # Parameter names used in the model functions are as in
+ # https://nbviewer.jupyter.org/urls/jrwb.de/nb/Symbolic%20ODE%20solutions%20for%20mkin.ipynb
+ types <- unname(sapply(mkin_model$spec, function(x) x$type))
+ if (length(mkin_model$spec) == 2 &! "SFORB" %in% types ) {
+ # Initial value for the metabolite (n20) must be fixed
+ if (names(odeini_fixed) == names(mkin_model$spec)[2]) {
+ n20 <- odeini_fixed
+ parent_name <- names(mkin_model$spec)[1]
+ if (identical(types, c("SFO", "SFO"))) {
+ if (transformations == "mkin") {
+ model_function <- function(psi, id, xidep) {
+ t <- xidep[, "time"]
+ n10 <- psi[id, 1]
+ k1 <- exp(psi[id, 2])
+ k2 <- exp(psi[id, 3])
+ f12 <- plogis(psi[id, 4])
+ ifelse(xidep[, "name"] == parent_name,
+ n10 * exp(- k1 * t),
+ (((k2 - k1) * n20 - f12 * k1 * n10) * exp(- k2 * t)) / (k2 - k1) +
+ (f12 * k1 * n10 * exp(- k1 * t)) / (k2 - k1)
+ )
+ }
+ } else {
+ model_function <- function(psi, id, xidep) {
+ t <- xidep[, "time"]
+ n10 <- psi[id, 1]
+ k1 <- psi[id, 2]
+ k2 <- psi[id, 3]
+ f12 <- psi[id, 4]
+ ifelse(xidep[, "name"] == parent_name,
+ n10 * exp(- k1 * t),
+ (((k2 - k1) * n20 - f12 * k1 * n10) * exp(- k2 * t)) / (k2 - k1) +
+ (f12 * k1 * n10 * exp(- k1 * t)) / (k2 - k1)
+ )
+ }
+ transform.par = c(0, 1, 1, 3)
+ }
+ }
+ if (identical(types, c("DFOP", "SFO"))) {
+ if (transformations == "mkin") {
+ model_function <- function(psi, id, xidep) {
+ t <- xidep[, "time"]
+ n10 <- psi[id, 1]
+ k2 <- exp(psi[id, 2])
+ f12 <- plogis(psi[id, 3])
+ l1 <- exp(psi[id, 4])
+ l2 <- exp(psi[id, 5])
+ g <- plogis(psi[id, 6])
+ ifelse(xidep[, "name"] == parent_name,
+ n10 * (g * exp(- l1 * t) + (1 - g) * exp(- l2 * t)),
+ ((f12 * g - f12) * l2 * n10 * exp(- l2 * t)) / (l2 - k2) -
+ (f12 * g * l1 * n10 * exp(- l1 * t)) / (l1 - k2) +
+ ((((l1 - k2) * l2 - k2 * l1 + k2^2) * n20 +
+ ((f12 * l1 + (f12 * g - f12) * k2) * l2 -
+ f12 * g * k2 * l1) * n10) * exp( - k2 * t)) /
+ ((l1 - k2) * l2 - k2 * l1 + k2^2)
+ )
+ }
+ } else {
+ model_function <- function(psi, id, xidep) {
+ t <- xidep[, "time"]
+ n10 <- psi[id, 1]
+ k2 <- psi[id, 2]
+ f12 <- psi[id, 3]
+ l1 <- psi[id, 4]
+ l2 <- psi[id, 5]
+ g <- psi[id, 6]
+ ifelse(xidep[, "name"] == parent_name,
+ n10 * (g * exp(- l1 * t) + (1 - g) * exp(- l2 * t)),
+ ((f12 * g - f12) * l2 * n10 * exp(- l2 * t)) / (l2 - k2) -
+ (f12 * g * l1 * n10 * exp(- l1 * t)) / (l1 - k2) +
+ ((((l1 - k2) * l2 - k2 * l1 + k2^2) * n20 +
+ ((f12 * l1 + (f12 * g - f12) * k2) * l2 -
+ f12 * g * k2 * l1) * n10) * exp( - k2 * t)) /
+ ((l1 - k2) * l2 - k2 * l1 + k2^2)
+ )
+ }
+ transform.par = c(0, 1, 3, 1, 1, 3)
+ }
+ }
+ }
+ }
+ }
+
+ if (is.function(model_function) & solution_type == "auto") {
+ solution_type = "analytical saemix"
+ } else {
+
+ if (solution_type == "auto")
+ solution_type <- object[[1]]$solution_type
+
+ model_function <- function(psi, id, xidep) {
+
+ uid <- unique(id)
+
+ res_list <- lapply(uid, function(i) {
+
+ transparms_optim <- as.numeric(psi[i, ]) # psi[i, ] is a dataframe when called in saemix.predict
+ names(transparms_optim) <- names(degparms_optim)
+
+ odeini_optim <- transparms_optim[odeini_optim_parm_names]
+ names(odeini_optim) <- gsub('_0$', '', odeini_optim_parm_names)
+
+ odeini <- c(odeini_optim, odeini_fixed)[names(mkin_model$diffs)]
+
+ ode_transparms_optim_names <- setdiff(names(transparms_optim), odeini_optim_parm_names)
+ odeparms_optim <- backtransform_odeparms(transparms_optim[ode_transparms_optim_names], mkin_model,
+ transform_rates = object[[1]]$transform_rates,
+ transform_fractions = object[[1]]$transform_fractions)
+ odeparms <- c(odeparms_optim, odeparms_fixed)
+
+ xidep_i <- subset(xidep, id == i)
+
+ if (solution_type == "analytical") {
+ out_values <- mkin_model$deg_func(xidep_i, odeini, odeparms)
+ } else {
+
+ i_time <- xidep_i$time
+ i_name <- xidep_i$name
+
+ out_wide <- mkinpredict(mkin_model,
+ odeparms = odeparms, odeini = odeini,
+ solution_type = solution_type,
+ outtimes = sort(unique(i_time)),
+ na_stop = FALSE
+ )
+
+ out_index <- cbind(as.character(i_time), as.character(i_name))
+ out_values <- out_wide[out_index]
+ }
+ return(out_values)
+ })
+ res <- unlist(res_list)
+ return(res)
+ }
+ }
+
+ error.model <- switch(object[[1]]$err_mod,
+ const = "constant",
+ tc = "combined",
+ obs = "constant")
+
+ if (object[[1]]$err_mod == "obs") {
+ warning("The error model 'obs' (variance by variable) can currently not be transferred to an saemix model")
+ }
+
+ error.init <- switch(object[[1]]$err_mod,
+ const = c(a = mean(sapply(object, function(x) x$errparms)), b = 1),
+ tc = c(a = mean(sapply(object, function(x) x$errparms[1])),
+ b = mean(sapply(object, function(x) x$errparms[2]))),
+ obs = c(a = mean(sapply(object, function(x) x$errparms)), b = 1))
+
+ degparms_psi0 <- degparms_optim
+ degparms_psi0[names(degparms_start)] <- degparms_start
+ psi0_matrix <- matrix(degparms_psi0, nrow = 1)
+ colnames(psi0_matrix) <- names(degparms_psi0)
+
+ res <- saemix::saemixModel(model_function,
+ psi0 = psi0_matrix,
+ "Mixed model generated from mmkin object",
+ transform.par = transform.par,
+ error.model = error.model,
+ verbose = verbose
+ )
+ attr(res, "mean_dp_start") <- degparms_optim
+ return(res)
+}
+
+#' @rdname saem
+#' @return An [saemix::SaemixData] object.
+#' @export
+saemix_data <- function(object, verbose = FALSE, ...) {
+ if (nrow(object) > 1) stop("Only row objects allowed")
+ ds_names <- colnames(object)
+
+ ds_list <- lapply(object, function(x) x$data[c("time", "variable", "observed")])
+ names(ds_list) <- ds_names
+ ds_saemix_all <- purrr::map_dfr(ds_list, function(x) x, .id = "ds")
+ ds_saemix <- data.frame(ds = ds_saemix_all$ds,
+ name = as.character(ds_saemix_all$variable),
+ time = ds_saemix_all$time,
+ value = ds_saemix_all$observed,
+ stringsAsFactors = FALSE)
+
+ res <- saemix::saemixData(ds_saemix,
+ name.group = "ds",
+ name.predictors = c("time", "name"),
+ name.response = "value",
+ verbose = verbose,
+ ...)
+ return(res)
+}
diff --git a/R/summary.saem.mmkin.R b/R/summary.saem.mmkin.R
new file mode 100644
index 00000000..e92c561c
--- /dev/null
+++ b/R/summary.saem.mmkin.R
@@ -0,0 +1,268 @@
+#' Summary method for class "saem.mmkin"
+#'
+#' Lists model equations, initial parameter values, optimised parameters
+#' for fixed effects (population), random effects (deviations from the
+#' population mean) and residual error model, as well as the resulting
+#' endpoints such as formation fractions and DT50 values. Optionally
+#' (default is FALSE), the data are listed in full.
+#'
+#' @param object an object of class [saem.mmkin]
+#' @param x an object of class [summary.saem.mmkin]
+#' @param data logical, indicating whether the full data should be included in
+#' the summary.
+#' @param verbose Should the summary be verbose?
+#' @param distimes logical, indicating whether DT50 and DT90 values should be
+#' included.
+#' @param digits Number of digits to use for printing
+#' @param \dots optional arguments passed to methods like \code{print}.
+#' @return The summary function returns a list based on the [saemix::SaemixObject]
+#' obtained in the fit, with at least the following additional components
+#' \item{saemixversion, mkinversion, Rversion}{The saemix, mkin and R versions used}
+#' \item{date.fit, date.summary}{The dates where the fit and the summary were
+#' produced}
+#' \item{diffs}{The differential equations used in the degradation model}
+#' \item{use_of_ff}{Was maximum or minimum use made of formation fractions}
+#' \item{data}{The data}
+#' \item{confint_trans}{Transformed parameters as used in the optimisation, with confidence intervals}
+#' \item{confint_back}{Backtransformed parameters, with confidence intervals if available}
+#' \item{confint_errmod}{Error model parameters with confidence intervals}
+#' \item{ff}{The estimated formation fractions derived from the fitted
+#' model.}
+#' \item{distimes}{The DT50 and DT90 values for each observed variable.}
+#' \item{SFORB}{If applicable, eigenvalues of SFORB components of the model.}
+#' The print method is called for its side effect, i.e. printing the summary.
+#' @importFrom stats predict vcov
+#' @author Johannes Ranke for the mkin specific parts
+#' saemix authors for the parts inherited from saemix.
+#' @examples
+#' # Generate five datasets following DFOP-SFO kinetics
+#' sampling_times = c(0, 1, 3, 7, 14, 28, 60, 90, 120)
+#' dfop_sfo <- mkinmod(parent = mkinsub("DFOP", "m1"),
+#' m1 = mkinsub("SFO"), quiet = TRUE)
+#' set.seed(1234)
+#' k1_in <- rlnorm(5, log(0.1), 0.3)
+#' k2_in <- rlnorm(5, log(0.02), 0.3)
+#' g_in <- plogis(rnorm(5, qlogis(0.5), 0.3))
+#' f_parent_to_m1_in <- plogis(rnorm(5, qlogis(0.3), 0.3))
+#' k_m1_in <- rlnorm(5, log(0.02), 0.3)
+#'
+#' pred_dfop_sfo <- function(k1, k2, g, f_parent_to_m1, k_m1) {
+#' mkinpredict(dfop_sfo,
+#' c(k1 = k1, k2 = k2, g = g, f_parent_to_m1 = f_parent_to_m1, k_m1 = k_m1),
+#' c(parent = 100, m1 = 0),
+#' sampling_times)
+#' }
+#'
+#' ds_mean_dfop_sfo <- lapply(1:5, function(i) {
+#' mkinpredict(dfop_sfo,
+#' c(k1 = k1_in[i], k2 = k2_in[i], g = g_in[i],
+#' f_parent_to_m1 = f_parent_to_m1_in[i], k_m1 = k_m1_in[i]),
+#' c(parent = 100, m1 = 0),
+#' sampling_times)
+#' })
+#' names(ds_mean_dfop_sfo) <- paste("ds", 1:5)
+#'
+#' ds_syn_dfop_sfo <- lapply(ds_mean_dfop_sfo, function(ds) {
+#' add_err(ds,
+#' sdfunc = function(value) sqrt(1^2 + value^2 * 0.07^2),
+#' n = 1)[[1]]
+#' })
+#'
+#' \dontrun{
+#' # Evaluate using mmkin and saem
+#' f_mmkin_dfop_sfo <- mmkin(list(dfop_sfo), ds_syn_dfop_sfo,
+#' quiet = TRUE, error_model = "tc", cores = 5)
+#' f_saem_dfop_sfo <- saem(f_mmkin_dfop_sfo)
+#' summary(f_saem_dfop_sfo, data = TRUE)
+#' }
+#'
+#' @export
+summary.saem.mmkin <- function(object, data = FALSE, verbose = FALSE, distimes = TRUE, ...) {
+
+ mod_vars <- names(object$mkinmod$diffs)
+
+ pnames <- names(object$mean_dp_start)
+ np <- length(pnames)
+
+ conf.int <- object$so@results@conf.int
+ rownames(conf.int) <- conf.int$name
+ confint_trans <- as.matrix(conf.int[pnames, c("estimate", "lower", "upper")])
+ colnames(confint_trans)[1] <- "est."
+
+ # In case objects were produced by earlier versions of saem
+ if (is.null(object$transformations)) object$transformations <- "mkin"
+
+ if (object$transformations == "mkin") {
+ bp <- backtransform_odeparms(confint_trans[, "est."], object$mkinmod,
+ object$transform_rates, object$transform_fractions)
+ bpnames <- names(bp)
+
+ # Transform boundaries of CI for one parameter at a time,
+ # with the exception of sets of formation fractions (single fractions are OK).
+ f_names_skip <- character(0)
+ for (box in mod_vars) { # Figure out sets of fractions to skip
+ f_names <- grep(paste("^f", box, sep = "_"), pnames, value = TRUE)
+ n_paths <- length(f_names)
+ if (n_paths > 1) f_names_skip <- c(f_names_skip, f_names)
+ }
+
+ confint_back <- matrix(NA, nrow = length(bp), ncol = 3,
+ dimnames = list(bpnames, colnames(confint_trans)))
+ confint_back[, "est."] <- bp
+
+ for (pname in pnames) {
+ if (!pname %in% f_names_skip) {
+ par.lower <- confint_trans[pname, "lower"]
+ par.upper <- confint_trans[pname, "upper"]
+ names(par.lower) <- names(par.upper) <- pname
+ bpl <- backtransform_odeparms(par.lower, object$mkinmod,
+ object$transform_rates,
+ object$transform_fractions)
+ bpu <- backtransform_odeparms(par.upper, object$mkinmod,
+ object$transform_rates,
+ object$transform_fractions)
+ confint_back[names(bpl), "lower"] <- bpl
+ confint_back[names(bpu), "upper"] <- bpu
+ }
+ }
+ } else {
+ confint_back <- confint_trans
+ }
+
+ # Correlation of fixed effects (inspired by summary.nlme)
+ varFix <- vcov(object$so)[1:np, 1:np]
+ stdFix <- sqrt(diag(varFix))
+ object$corFixed <- array(
+ t(varFix/stdFix)/stdFix,
+ dim(varFix),
+ list(pnames, pnames))
+
+ # Random effects
+ rnames <- paste0("SD.", pnames)
+ confint_ranef <- as.matrix(conf.int[rnames, c("estimate", "lower", "upper")])
+ colnames(confint_ranef)[1] <- "est."
+
+ # Error model
+ enames <- if (object$err_mod == "const") "a.1" else c("a.1", "b.1")
+ confint_errmod <- as.matrix(conf.int[enames, c("estimate", "lower", "upper")])
+ colnames(confint_errmod)[1] <- "est."
+
+
+ object$confint_trans <- confint_trans
+ object$confint_ranef <- confint_ranef
+ object$confint_errmod <- confint_errmod
+ object$confint_back <- confint_back
+
+ object$date.summary = date()
+ object$use_of_ff = object$mkinmod$use_of_ff
+ object$error_model_algorithm = object$mmkin_orig[[1]]$error_model_algorithm
+ err_mod = object$mmkin_orig[[1]]$err_mod
+
+ object$diffs <- object$mkinmod$diffs
+ object$print_data <- data # boolean: Should we print the data?
+ so_pred <- object$so@results@predictions
+
+ names(object$data)[4] <- "observed" # rename value to observed
+
+ object$verbose <- verbose
+
+ object$fixed <- object$mmkin_orig[[1]]$fixed
+ object$AIC = AIC(object$so)
+ object$BIC = BIC(object$so)
+ object$logLik = logLik(object$so, method = "is")
+
+ ep <- endpoints(object)
+ if (length(ep$ff) != 0)
+ object$ff <- ep$ff
+ if (distimes) object$distimes <- ep$distimes
+ if (length(ep$SFORB) != 0) object$SFORB <- ep$SFORB
+ class(object) <- c("summary.saem.mmkin")
+ return(object)
+}
+
+#' @rdname summary.saem.mmkin
+#' @export
+print.summary.saem.mmkin <- function(x, digits = max(3, getOption("digits") - 3), verbose = x$verbose, ...) {
+ cat("saemix version used for fitting: ", x$saemixversion, "\n")
+ cat("mkin version used for pre-fitting: ", x$mkinversion, "\n")
+ cat("R version used for fitting: ", x$Rversion, "\n")
+
+ cat("Date of fit: ", x$date.fit, "\n")
+ cat("Date of summary:", x$date.summary, "\n")
+
+ cat("\nEquations:\n")
+ nice_diffs <- gsub("^(d.*) =", "\\1/dt =", x[["diffs"]])
+ writeLines(strwrap(nice_diffs, exdent = 11))
+
+ cat("\nData:\n")
+ cat(nrow(x$data), "observations of",
+ length(unique(x$data$name)), "variable(s) grouped in",
+ length(unique(x$data$ds)), "datasets\n")
+
+ cat("\nModel predictions using solution type", x$solution_type, "\n")
+
+ cat("\nFitted in", x$time[["elapsed"]], "s using", paste(x$so@options$nbiter.saemix, collapse = ", "), "iterations\n")
+
+ cat("\nVariance model: ")
+ cat(switch(x$err_mod,
+ const = "Constant variance",
+ obs = "Variance unique to each observed variable",
+ tc = "Two-component variance function"), "\n")
+
+ cat("\nMean of starting values for individual parameters:\n")
+ print(x$mean_dp_start, digits = digits)
+
+ cat("\nFixed degradation parameter values:\n")
+ if(length(x$fixed$value) == 0) cat("None\n")
+ else print(x$fixed, digits = digits)
+
+ cat("\nResults:\n\n")
+ cat("Likelihood computed by importance sampling\n")
+ print(data.frame(AIC = x$AIC, BIC = x$BIC, logLik = x$logLik,
+ row.names = " "), digits = digits)
+
+ cat("\nOptimised parameters:\n")
+ print(x$confint_trans, digits = digits)
+
+ if (nrow(x$confint_trans) > 1) {
+ corr <- x$corFixed
+ class(corr) <- "correlation"
+ print(corr, title = "\nCorrelation:", ...)
+ }
+
+ cat("\nRandom effects:\n")
+ print(x$confint_ranef, digits = digits)
+
+ cat("\nVariance model:\n")
+ print(x$confint_errmod, digits = digits)
+
+ if (x$transformations == "mkin") {
+ cat("\nBacktransformed parameters:\n")
+ print(x$confint_back, digits = digits)
+ }
+
+ printSFORB <- !is.null(x$SFORB)
+ if(printSFORB){
+ cat("\nEstimated Eigenvalues of SFORB model(s):\n")
+ print(x$SFORB, digits = digits,...)
+ }
+
+ printff <- !is.null(x$ff)
+ if(printff){
+ cat("\nResulting formation fractions:\n")
+ print(data.frame(ff = x$ff), digits = digits,...)
+ }
+
+ printdistimes <- !is.null(x$distimes)
+ if(printdistimes){
+ cat("\nEstimated disappearance times:\n")
+ print(x$distimes, digits = digits,...)
+ }
+
+ if (x$print_data){
+ cat("\nData:\n")
+ print(format(x$data, digits = digits, ...), row.names = FALSE)
+ }
+
+ invisible(x)
+}
diff --git a/_pkgdown.yml b/_pkgdown.yml
index f632ddb0..5cfaeedf 100644
--- a/_pkgdown.yml
+++ b/_pkgdown.yml
@@ -1,8 +1,8 @@
url: https://pkgdown.jrwb.de/mkin
development:
- mode: release
- version_label: default
+ mode: devel
+ version_label: info
template:
bootswatch: spacelab
@@ -39,10 +39,13 @@ reference:
desc: Create and work with nonlinear mixed effects models
contents:
- nlme.mmkin
+ - saem.mmkin
- plot.mixed.mmkin
- summary.nlme.mmkin
+ - summary.saem.mmkin
- nlme_function
- get_deg_func
+ - saemix_model
- mixed
- title: Datasets and known results
contents:
diff --git a/build.log b/build.log
index 09d94134..7a3b024e 100644
--- a/build.log
+++ b/build.log
@@ -6,5 +6,5 @@
* creating vignettes ... OK
* checking for LF line-endings in source and make files and shell scripts
* checking for empty or unneeded directories
-* building ‘mkin_1.0.1.tar.gz’
+* building ‘mkin_1.0.1.9000.tar.gz’
diff --git a/check.log b/check.log
index d0aeb08a..58d4d6c7 100644
--- a/check.log
+++ b/check.log
@@ -5,10 +5,15 @@
* using options ‘--no-tests --as-cran’
* checking for file ‘mkin/DESCRIPTION’ ... OK
* checking extension type ... Package
-* this is package ‘mkin’ version ‘1.0.1’
+* this is package ‘mkin’ version ‘1.0.1.9000’
* package encoding: UTF-8
-* checking CRAN incoming feasibility ... Note_to_CRAN_maintainers
+* checking CRAN incoming feasibility ... NOTE
Maintainer: ‘Johannes Ranke ’
+
+Version contains large components (1.0.1.9000)
+
+Unknown, possibly mis-spelled, fields in DESCRIPTION:
+ ‘Remotes’
* checking package namespace information ... OK
* checking package dependencies ... OK
* checking if this is a source package ... OK
@@ -67,5 +72,9 @@ Maintainer: ‘Johannes Ranke ’
* checking for detritus in the temp directory ... OK
* DONE
-Status: OK
+Status: 1 NOTE
+See
+ ‘/home/jranke/git/mkin/mkin.Rcheck/00check.log’
+for details.
+
diff --git a/docs/dev/404.html b/docs/dev/404.html
index bea38406..5f29faf2 100644
--- a/docs/dev/404.html
+++ b/docs/dev/404.html
@@ -71,7 +71,7 @@
mkin
- 0.9.50.4
+ 1.0.1.9000
diff --git a/docs/dev/articles/index.html b/docs/dev/articles/index.html
index 6daa6960..441d49c0 100644
--- a/docs/dev/articles/index.html
+++ b/docs/dev/articles/index.html
@@ -71,7 +71,7 @@
mkin
- 0.9.50.4
+ 1.0.1.9000
diff --git a/docs/dev/authors.html b/docs/dev/authors.html
index d592b39f..9641eec0 100644
--- a/docs/dev/authors.html
+++ b/docs/dev/authors.html
@@ -71,7 +71,7 @@
mkin
- 0.9.50.4
+ 1.0.1.9000
@@ -147,15 +147,15 @@
Katrin Lindenberger. Contributor.
-
+
contributed to mkinresplot()
René Lehmann. Contributor.
-
+
ilr() and invilr()
Eurofins Regulatory AG. Copyright holder.
-
+
copyright for some of the contributions of JR 2012-2014
diff --git a/docs/dev/index.html b/docs/dev/index.html
index a4399963..8888633d 100644
--- a/docs/dev/index.html
+++ b/docs/dev/index.html
@@ -38,7 +38,7 @@
mkin
- 0.9.50.4
+ 1.0.1.9000
@@ -147,7 +147,7 @@
Three different error models can be selected using the argument error_model
to the mkinfit
function.
The ‘variance by variable’ error model which is often fitted using Iteratively Reweighted Least Squares (IRLS) should now be specified as error_model = "obs"
.
A two-component error model similar to the one proposed by Rocke and Lorenzato can be selected using the argument error_model = "tc"
.
-Nonlinear mixed-effects models can be created from fits of the same degradation model to different datasets for the same compound by using the nlme.mmkin method.
+Nonlinear mixed-effects models can be created from fits of the same degradation model to different datasets for the same compound by using the nlme.mmkin method. Note that the convergence of the nlme fits depends on the quality of the data. Convergence is better for simple models and data for many groups (e.g. soils).
diff --git a/docs/dev/news/index.html b/docs/dev/news/index.html
index c3597efe..998917f2 100644
--- a/docs/dev/news/index.html
+++ b/docs/dev/news/index.html
@@ -71,7 +71,7 @@
mkin
- 0.9.50.4
+ 1.0.1.9000
@@ -141,13 +141,36 @@
Source: NEWS.md
-
-