diff options
Diffstat (limited to 'R/saemix.R')
-rw-r--r-- | R/saemix.R | 42 |
1 files changed, 41 insertions, 1 deletions
@@ -25,7 +25,7 @@ #' \dontrun{ #' f_mmkin <- mmkin(list("SFO-SFO" = sfo_sfo), ds, quiet = TRUE) #' library(saemix) -#' m_saemix <- saemix_model(f_mmkin) +#' m_saemix <- saemix_model(f_mmkin, cores = 1) #' d_saemix <- saemix_data(f_mmkin) #' saemix_options <- list(seed = 123456, #' save = FALSE, save.graphs = FALSE, displayProgress = FALSE, @@ -33,6 +33,30 @@ #' f_saemix <- saemix(m_saemix, d_saemix, saemix_options) #' plot(f_saemix, plot.type = "convergence") #' } +#' # Synthetic data with two-component error +#' sampling_times = c(0, 1, 3, 7, 14, 28, 60, 90, 120) +#' dt50_sfo_in <- c(800, 900, 1000, 1111.11, 1250) +#' k_in <- log(2) / dt50_sfo_in +#' +#' SFO <- mkinmod(parent = mkinsub("SFO")) +#' +#' pred_sfo <- function(k) { +#' mkinpredict(SFO, c(k_parent = k), +#' c(parent = 100), sampling_times) +#' } +#' +#' ds_sfo_mean <- lapply(k_in, pred_sfo) +#' set.seed(123456L) +#' ds_sfo_syn <- lapply(ds_sfo_mean, function(ds) { +#' add_err(ds, sdfunc = function(value) sqrt(1^2 + value^2 * 0.07^2), +#' n = 1)[[1]] +#' }) +#' f_mmkin_syn <- mmkin("SFO", ds_sfo_syn, error_model = "tc", quiet = TRUE) +#' m_saemix_tc <- saemix_model(f_mmkin_syn, cores = 1) +#' d_saemix_tc <- saemix_data(f_mmkin_syn) +#' f_saemix_tc <- saemix(m_saemix_tc, d_saemix_tc, saemix_options) +#' plot(f_saemix_tc, plot.type = "convergence") +#' #' @return An [saemix::SaemixModel] object. #' @export saemix_model <- function(object, cores = parallel::detectCores()) { @@ -101,9 +125,25 @@ saemix_model <- function(object, cores = parallel::detectCores()) { raneff_0 <- mean_degparms(object, random = TRUE)$random$ds var_raneff_0 <- apply(raneff_0, 2, var) + 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) was not transferred to the 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)) + res <- saemixModel(model_function, psi0, "Mixed model generated from mmkin object", transform.par = rep(0, length(degparms_optim)), + error.model = error.model, error.init = error.init, omega.init = diag(var_raneff_0) ) return(res) |