From 272aba066f0d7502e319b7e7f14009318cd44348 Mon Sep 17 00:00:00 2001 From: Johannes Ranke Date: Thu, 15 Oct 2020 12:27:40 +0200 Subject: Updating docs of current dev version --- R/saemix.R | 42 ++- docs/dev/404.html | 4 +- docs/dev/articles/index.html | 4 +- docs/dev/articles/mkin.html | 45 +-- .../mkin_files/figure-html/unnamed-chunk-2-1.png | Bin 116140 -> 92716 bytes docs/dev/articles/web_only/NAFTA_examples.html | 154 ++++++---- .../NAFTA_examples_files/figure-html/p10-1.png | Bin 53828 -> 80950 bytes .../NAFTA_examples_files/figure-html/p11-1.png | Bin 51210 -> 77737 bytes .../NAFTA_examples_files/figure-html/p12a-1.png | Bin 54489 -> 83111 bytes .../NAFTA_examples_files/figure-html/p12b-1.png | Bin 47070 -> 72088 bytes .../NAFTA_examples_files/figure-html/p13-1.png | Bin 51343 -> 79358 bytes .../NAFTA_examples_files/figure-html/p14-1.png | Bin 53898 -> 82078 bytes .../NAFTA_examples_files/figure-html/p15a-1.png | Bin 51113 -> 78259 bytes .../NAFTA_examples_files/figure-html/p15b-1.png | Bin 52979 -> 80090 bytes .../NAFTA_examples_files/figure-html/p16-1.png | Bin 62414 -> 95465 bytes .../NAFTA_examples_files/figure-html/p5a-1.png | Bin 55286 -> 84037 bytes .../NAFTA_examples_files/figure-html/p5b-1.png | Bin 53505 -> 81732 bytes .../NAFTA_examples_files/figure-html/p6-1.png | Bin 54907 -> 84055 bytes .../NAFTA_examples_files/figure-html/p7-1.png | Bin 66709 -> 103768 bytes .../NAFTA_examples_files/figure-html/p8-1.png | Bin 61400 -> 94196 bytes .../NAFTA_examples_files/figure-html/p9a-1.png | Bin 53005 -> 80132 bytes .../NAFTA_examples_files/figure-html/p9b-1.png | Bin 49914 -> 77691 bytes docs/dev/authors.html | 4 +- docs/dev/index.html | 7 +- docs/dev/news/index.html | 16 +- docs/dev/pkgdown.css | 4 +- docs/dev/pkgdown.yml | 12 +- docs/dev/reference/endpoints.html | 30 +- docs/dev/reference/index.html | 4 +- docs/dev/reference/saemix-1.png | Bin 31551 -> 50490 bytes docs/dev/reference/saemix-2.png | Bin 58815 -> 60503 bytes docs/dev/reference/saemix-3.png | Bin 0 -> 16942 bytes docs/dev/reference/saemix-4.png | Bin 0 -> 65995 bytes docs/dev/reference/saemix.html | 323 ++++++++++++++++++--- docs/dev/sitemap.xml | 3 + man/saemix.Rd | 26 +- 36 files changed, 522 insertions(+), 156 deletions(-) create mode 100644 docs/dev/reference/saemix-3.png create mode 100644 docs/dev/reference/saemix-4.png diff --git a/R/saemix.R b/R/saemix.R index 24c0f0d0..eb64bf79 100644 --- a/R/saemix.R +++ b/R/saemix.R @@ -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) diff --git a/docs/dev/404.html b/docs/dev/404.html index 9d14af1e..8ca69469 100644 --- a/docs/dev/404.html +++ b/docs/dev/404.html @@ -71,7 +71,7 @@ mkin - 0.9.50.3 + 0.9.50.4 @@ -160,7 +160,7 @@ Content not found. Please use links in the navbar.
-

Site built with pkgdown 1.5.1.

+

Site built with pkgdown 1.6.1.

diff --git a/docs/dev/articles/index.html b/docs/dev/articles/index.html index b8ba57d8..81d49aa3 100644 --- a/docs/dev/articles/index.html +++ b/docs/dev/articles/index.html @@ -71,7 +71,7 @@ mkin - 0.9.50.3 + 0.9.50.4 @@ -173,7 +173,7 @@
-

Site built with pkgdown 1.5.1.

+

Site built with pkgdown 1.6.1.

diff --git a/docs/dev/articles/mkin.html b/docs/dev/articles/mkin.html index 6865fe96..07627029 100644 --- a/docs/dev/articles/mkin.html +++ b/docs/dev/articles/mkin.html @@ -32,7 +32,7 @@ mkin - 0.9.50.3 + 0.9.50.4 @@ -80,7 +80,7 @@
  • - + @@ -146,7 +146,7 @@ advantage that the SFORB model can also be used for metabolites." />

    Function to calculate endpoints for further use from kinetic models fitted with mkinfit

    - Source: R/endpoints.R + Source: R/endpoints.R
    @@ -158,7 +158,7 @@ are equivalent to the rate constants of the DFOP model, but with the advantage that the SFORB model can also be used for metabolites.

    - 
    endpoints(fit)
    + 
    endpoints(fit)

    Arguments

    @@ -179,20 +179,26 @@ of these SFORB models, equivalent to DFOP rate constants

    Note

    The function is used internally by summary.mkinfit.

    +

    Author

    + +

    Johannes Ranke

    Examples

    
-  fit <- mkinfit("FOMC", FOCUS_2006_C, quiet = TRUE)
-  endpoints(fit)
    #> $distimes
+  fit <- mkinfit("FOMC", FOCUS_2006_C, quiet = TRUE)
+  endpoints(fit)
+
    #> $distimes
 #>            DT50    DT90 DT50back
 #> parent 1.785233 15.1479 4.559973
 #> 
      # \dontrun{
-    fit_2 <- mkinfit("DFOP", FOCUS_2006_C, quiet = TRUE)
-    endpoints(fit_2)
    #> $distimes
+    fit_2 <- mkinfit("DFOP", FOCUS_2006_C, quiet = TRUE)
+    endpoints(fit_2)
+
    #> $distimes
 #>            DT50     DT90 DT50back  DT50_k1  DT50_k2
 #> parent 1.886925 21.25106 6.397207 1.508293 38.83438
-#> 
        fit_3 <- mkinfit("SFORB", FOCUS_2006_C, quiet = TRUE)
-    endpoints(fit_3)
    #> $ff
+#> 
        fit_3 <- mkinfit("SFORB", FOCUS_2006_C, quiet = TRUE)
+    endpoints(fit_3)
+
    #> $ff
 #> parent_free 
 #>           1 
 #> 
@@ -203,7 +209,7 @@ of these SFORB models, equivalent to DFOP rate constants
    
 #> $distimes
 #>            DT50     DT90 DT50back DT50_parent_b1 DT50_parent_b2
 #> parent 1.886925 21.25106 6.397208       1.508293       38.83438
-#> 
      # }
+#> 
      # }
    @@ -221,7 +227,7 @@ of these SFORB models, equivalent to DFOP rate constants

    -

    Site built with pkgdown 1.5.1.

    +

    Site built with pkgdown 1.6.1.

    diff --git a/docs/dev/reference/index.html b/docs/dev/reference/index.html index a5cd6663..5d70a450 100644 --- a/docs/dev/reference/index.html +++ b/docs/dev/reference/index.html @@ -71,7 +71,7 @@ mkin - 0.9.50.3 + 0.9.50.4 @@ -678,7 +678,7 @@ kinetic models fitted with mkinfit

    -

    Site built with pkgdown 1.5.1.

    +

    Site built with pkgdown 1.6.1.

    diff --git a/docs/dev/reference/saemix-1.png b/docs/dev/reference/saemix-1.png index 0d79300d..696fa7d9 100644 Binary files a/docs/dev/reference/saemix-1.png and b/docs/dev/reference/saemix-1.png differ diff --git a/docs/dev/reference/saemix-2.png b/docs/dev/reference/saemix-2.png index 04de70b5..6525a36c 100644 Binary files a/docs/dev/reference/saemix-2.png and b/docs/dev/reference/saemix-2.png differ diff --git a/docs/dev/reference/saemix-3.png b/docs/dev/reference/saemix-3.png new file mode 100644 index 00000000..13157f91 Binary files /dev/null and b/docs/dev/reference/saemix-3.png differ diff --git a/docs/dev/reference/saemix-4.png b/docs/dev/reference/saemix-4.png new file mode 100644 index 00000000..1e717cb2 Binary files /dev/null and b/docs/dev/reference/saemix-4.png differ diff --git a/docs/dev/reference/saemix.html b/docs/dev/reference/saemix.html index da28bff1..5b160959 100644 --- a/docs/dev/reference/saemix.html +++ b/docs/dev/reference/saemix.html @@ -75,7 +75,7 @@ a list of datasets." /> mkin - 0.9.50.3 + 0.9.50.4 @@ -123,7 +123,7 @@ a list of datasets." />
    • - + @@ -142,7 +142,7 @@ a list of datasets." />

      Create saemix models from mmkin row objects

      - Source: R/saemix.R + Source: R/saemix.R
      @@ -153,9 +153,9 @@ list of mkinfit objects that have been obtained by fitting the same model to a list of datasets.

      - 
      saemix_model(object, cores = parallel::detectCores())
+    
      saemix_model(object, cores = parallel::detectCores())
 
-saemix_data(object, ...)
      
+saemix_data(object, ...)

      Arguments

    @@ -187,14 +187,17 @@ mmkin. Starting variances of the random effects (argument omega.init) are the variances of the deviations of the parameters from these mean values.

    Examples

    - 
    ds <- lapply(experimental_data_for_UBA_2019[6:10],
- function(x) subset(x$data[c("name", "time", "value")]))
-names(ds) <- paste("Dataset", 6:10)
-sfo_sfo <- mkinmod(parent = mkinsub("SFO", "A1"),
-  A1 = mkinsub("SFO"))
    #> Successfully compiled differential equation model from auto-generated C code.
    # \dontrun{
-f_mmkin <- mmkin(list("SFO-SFO" = sfo_sfo), ds, quiet = TRUE)
-library(saemix)
    #> Package saemix, version 3.1.9000
-#>   please direct bugs, questions and feedback to emmanuelle.comets@inserm.fr
    m_saemix <- saemix_model(f_mmkin)
    #> 
+    
    ds <- lapply(experimental_data_for_UBA_2019[6:10],
+ function(x) subset(x$data[c("name", "time", "value")]))
+names(ds) <- paste("Dataset", 6:10)
+sfo_sfo <- mkinmod(parent = mkinsub("SFO", "A1"),
+  A1 = mkinsub("SFO"))
+
    #> Successfully compiled differential equation model from auto-generated C code.
    # \dontrun{
+f_mmkin <- mmkin(list("SFO-SFO" = sfo_sfo), ds, quiet = TRUE)
+library(saemix)
+
    #> Package saemix, version 3.1.9000
+#>   please direct bugs, questions and feedback to emmanuelle.comets@inserm.fr
    m_saemix <- saemix_model(f_mmkin, cores = 1)
+
    #> 
 #> 
 #> The following SaemixModel object was successfully created:
 #> 
@@ -234,8 +237,8 @@ variances of the deviations of the parameters from these mean values.
    
 #>     res <- unlist(res_list)
 #>     return(res)
 #> }
-#> <bytecode: 0x55555bbf17c0>
-#> <environment: 0x55555bbec4e0>
+#> <bytecode: 0x55555c1e7720>
+#> <environment: 0x55555c1eff38>
 #>   Nb of parameters: 4 
 #>       parameter names:  parent_0 log_k_parent log_k_A1 f_parent_ilr_1 
 #>       distribution:
@@ -250,11 +253,12 @@ variances of the deviations of the parameters from these mean values.
    
 #> log_k_parent          0            1        0              0
 #> log_k_A1              0            0        1              0
 #> f_parent_ilr_1        0            0        0              1
-#>   Error model: constant , initial values: a.1=1 
+#>   Error model: constant , initial values: a.1=4.97259024646577 
 #>     No covariate in the model.
 #>     Initial values
 #>              parent_0 log_k_parent  log_k_A1 f_parent_ilr_1
-#> Pop.CondInit 86.53449    -3.207005 -3.060308      -1.920449
    d_saemix <- saemix_data(f_mmkin)
    #> 
+#> Pop.CondInit 86.53449    -3.207005 -3.060308      -1.920449
    d_saemix <- saemix_data(f_mmkin)
+
    #> 
 #> 
 #> The following SaemixData object was successfully created:
 #> 
@@ -262,14 +266,15 @@ variances of the deviations of the parameters from these mean values.
    
 #>     longitudinal data for use with the SAEM algorithm
 #> Dataset ds_saemix 
 #>     Structured data: value ~ time + name | ds 
-#>     X variable for graphs: time () 
    saemix_options <- list(seed = 123456,
-  save = FALSE, save.graphs = FALSE, displayProgress = FALSE,
-  nbiter.saemix = c(200, 80))
-f_saemix <- saemix(m_saemix, d_saemix, saemix_options)
    #> Running main SAEM algorithm
-#> [1] "Wed May 27 07:45:07 2020"
+#>     X variable for graphs: time () 
    saemix_options <- list(seed = 123456,
+  save = FALSE, save.graphs = FALSE, displayProgress = FALSE,
+  nbiter.saemix = c(200, 80))
+f_saemix <- saemix(m_saemix, d_saemix, saemix_options)
+
    #> Running main SAEM algorithm
+#> [1] "Thu Oct 15 11:55:14 2020"
 #> ..
 #>     Minimisation finished
-#> [1] "Wed May 27 07:51:24 2020"
    #> Nonlinear mixed-effects model fit by the SAEM algorithm
+#> [1] "Thu Oct 15 11:57:07 2020"
    #> Nonlinear mixed-effects model fit by the SAEM algorithm
 #> -----------------------------------
 #> ----          Data             ----
 #> -----------------------------------
@@ -333,8 +338,8 @@ variances of the deviations of the parameters from these mean values.
    
 #>     res <- unlist(res_list)
 #>     return(res)
 #> }
-#> <bytecode: 0x55555bbf17c0>
-#> <environment: 0x55555bbec4e0>
+#> <bytecode: 0x55555c1e7720>
+#> <environment: 0x55555c1eff38>
 #>   Nb of parameters: 4 
 #>       parameter names:  parent_0 log_k_parent log_k_A1 f_parent_ilr_1 
 #>       distribution:
@@ -349,7 +354,7 @@ variances of the deviations of the parameters from these mean values.
    
 #> log_k_parent          0            1        0              0
 #> log_k_A1              0            0        1              0
 #> f_parent_ilr_1        0            0        0              1
-#>   Error model: constant , initial values: a.1=1 
+#>   Error model: constant , initial values: a.1=4.97259024646577 
 #>     No covariate in the model.
 #>     Initial values
 #>              parent_0 log_k_parent  log_k_A1 f_parent_ilr_1
@@ -375,20 +380,20 @@ variances of the deviations of the parameters from these mean values.
    
 #> ----------------------------------------------------
 #> -----------------  Fixed effects  ------------------
 #> ----------------------------------------------------
-#>      Parameter      Estimate SE   CV(%)
-#> [1,] parent_0       86.14    1.61  1.9 
-#> [2,] log_k_parent   -3.21    0.59 18.5 
-#> [3,] log_k_A1       -4.66    0.30  6.4 
-#> [4,] f_parent_ilr_1 -0.33    0.30 91.7 
-#> [5,] a.1             4.68    0.27  5.8 
+#>   Parameter      Estimate SE   CV(%)
+#>   parent_0       86.09    1.57  1.8 
+#>   log_k_parent   -3.21    0.59 18.5 
+#>   log_k_A1       -4.69    0.31  6.6 
+#>   f_parent_ilr_1 -0.34    0.30 89.2 
+#> a a.1             4.69    0.27  5.8 
 #> ----------------------------------------------------
 #> -----------  Variance of random effects  -----------
 #> ----------------------------------------------------
 #>                Parameter             Estimate SE   CV(%)
-#> parent_0       omega2.parent_0       7.71     8.14 106  
-#> log_k_parent   omega2.log_k_parent   1.76     1.12  63  
-#> log_k_A1       omega2.log_k_A1       0.26     0.26 101  
-#> f_parent_ilr_1 omega2.f_parent_ilr_1 0.39     0.28  71  
+#> parent_0       omega2.parent_0       7.07     7.72 109  
+#> log_k_parent   omega2.log_k_parent   1.75     1.11  63  
+#> log_k_A1       omega2.log_k_A1       0.28     0.28  99  
+#> f_parent_ilr_1 omega2.f_parent_ilr_1 0.39     0.27  71  
 #> ----------------------------------------------------
 #> ------  Correlation matrix of random effects  ------
 #> ----------------------------------------------------
@@ -406,15 +411,245 @@ variances of the deviations of the parameters from these mean values.
    
 #> ---------------  Statistical criteria  -------------
 #> ----------------------------------------------------
 #> Likelihood computed by linearisation
-#>       -2LL= 1064.364 
-#>       AIC = 1082.364 
-#>       BIC = 1078.848 
+#>       -2LL= 1064.35 
+#>       AIC = 1082.35 
+#>       BIC = 1078.835 
 #> 
 #> Likelihood computed by importance sampling
-#>       -2LL= 1063.462 
-#>       AIC = 1081.462 
-#>       BIC = 1077.947 
-#> ----------------------------------------------------
    plot(f_saemix, plot.type = "convergence")
    #> Plotting convergence plots
    # }
+#>       -2LL= 1063.475 
+#>       AIC = 1081.475 
+#>       BIC = 1077.96 
+#> ----------------------------------------------------
    plot(f_saemix, plot.type = "convergence")
+
    #> Plotting convergence plots
    # }
+# 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)
+
    #> 
+#> 
+#> The following SaemixModel object was successfully created:
+#> 
+#> Nonlinear mixed-effects model
+#>   Model function:  Mixed model generated from mmkin object  Model type:  structural
+#> function (psi, id, xidep) 
+#> {
+#>     uid <- unique(id)
+#>     res_list <- parallel::mclapply(uid, function(i) {
+#>         transparms_optim <- psi[i, ]
+#>         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 (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 = object[[1]]$solution_type, 
+#>                 outtimes = sort(unique(i_time)))
+#>             out_index <- cbind(as.character(i_time), as.character(i_name))
+#>             out_values <- out_wide[out_index]
+#>         }
+#>         return(out_values)
+#>     }, mc.cores = cores)
+#>     res <- unlist(res_list)
+#>     return(res)
+#> }
+#> <bytecode: 0x55555c1e7720>
+#> <environment: 0x55555be5b848>
+#>   Nb of parameters: 2 
+#>       parameter names:  parent_0 log_k_parent 
+#>       distribution:
+#>      Parameter    Distribution Estimated
+#> [1,] parent_0     normal       Estimated
+#> [2,] log_k_parent normal       Estimated
+#>   Variance-covariance matrix:
+#>              parent_0 log_k_parent
+#> parent_0            1            0
+#> log_k_parent        0            1
+#>   Error model: combined , initial values: a.1=2.4206146215511 b.1=0.03831266409084 
+#>     No covariate in the model.
+#>     Initial values
+#>              parent_0 log_k_parent
+#> Pop.CondInit 100.2498     -9.33922
    d_saemix_tc <- saemix_data(f_mmkin_syn)
+
    #> 
+#> 
+#> The following SaemixData object was successfully created:
+#> 
+#> Object of class SaemixData
+#>     longitudinal data for use with the SAEM algorithm
+#> Dataset ds_saemix 
+#>     Structured data: value ~ time + name | ds 
+#>     X variable for graphs: time () 
    f_saemix_tc <- saemix(m_saemix_tc, d_saemix_tc, saemix_options)
+
    #> Running main SAEM algorithm
+#> [1] "Thu Oct 15 11:57:25 2020"
+#> ..
+#>     Minimisation finished
+#> [1] "Thu Oct 15 11:58:16 2020"
+#> Error in solve.default(Fth) : 
+#>   Lapack routine dgesv: system is exactly singular: U[2,2] = 0
+#> Error computing the Fisher Information Matrix: singular system.
+#> Error in solve.default(FO) : 
+#>   Lapack routine dgesv: system is exactly singular: U[2,2] = 0
+#> Error computing the Fisher Information Matrix: singular system.
    #> Nonlinear mixed-effects model fit by the SAEM algorithm
+#> -----------------------------------
+#> ----          Data             ----
+#> -----------------------------------
+#> Object of class SaemixData
+#>     longitudinal data for use with the SAEM algorithm
+#> Dataset ds_saemix 
+#>     Structured data: value ~ time + name | ds 
+#>     X variable for graphs: time () 
+#> Dataset characteristics:
+#>     number of subjects:     5 
+#>     number of observations: 90 
+#>     average/min/max nb obs: 18.00  /  18  /  18 
+#> First 10 lines of data:
+#>    ds time   name value mdv cens occ ytype
+#> 1   1    0 parent 105.9   0    0   1     1
+#> 2   1    0 parent  98.0   0    0   1     1
+#> 3   1    1 parent  97.4   0    0   1     1
+#> 4   1    1 parent 100.5   0    0   1     1
+#> 5   1    3 parent 115.6   0    0   1     1
+#> 6   1    3 parent 105.6   0    0   1     1
+#> 7   1    7 parent 108.6   0    0   1     1
+#> 8   1    7 parent 117.0   0    0   1     1
+#> 9   1   14 parent 107.0   0    0   1     1
+#> 10  1   14 parent  95.8   0    0   1     1
+#> -----------------------------------
+#> ----          Model            ----
+#> -----------------------------------
+#> Nonlinear mixed-effects model
+#>   Model function:  Mixed model generated from mmkin object  Model type:  structural
+#> function (psi, id, xidep) 
+#> {
+#>     uid <- unique(id)
+#>     res_list <- parallel::mclapply(uid, function(i) {
+#>         transparms_optim <- psi[i, ]
+#>         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 (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 = object[[1]]$solution_type, 
+#>                 outtimes = sort(unique(i_time)))
+#>             out_index <- cbind(as.character(i_time), as.character(i_name))
+#>             out_values <- out_wide[out_index]
+#>         }
+#>         return(out_values)
+#>     }, mc.cores = cores)
+#>     res <- unlist(res_list)
+#>     return(res)
+#> }
+#> <bytecode: 0x55555c1e7720>
+#> <environment: 0x55555be5b848>
+#>   Nb of parameters: 2 
+#>       parameter names:  parent_0 log_k_parent 
+#>       distribution:
+#>      Parameter    Distribution Estimated
+#> [1,] parent_0     normal       Estimated
+#> [2,] log_k_parent normal       Estimated
+#>   Variance-covariance matrix:
+#>              parent_0 log_k_parent
+#> parent_0            1            0
+#> log_k_parent        0            1
+#>   Error model: combined , initial values: a.1=2.4206146215511 b.1=0.03831266409084 
+#>     No covariate in the model.
+#>     Initial values
+#>              parent_0 log_k_parent
+#> Pop.CondInit 100.2498     -9.33922
+#> -----------------------------------
+#> ----    Key algorithm options  ----
+#> -----------------------------------
+#>     Estimation of individual parameters (MAP)
+#>     Estimation of standard errors and linearised log-likelihood
+#>     Estimation of log-likelihood by importance sampling
+#>     Number of iterations:  K1=200, K2=80 
+#>     Number of chains:  10 
+#>     Seed:  123456 
+#>     Number of MCMC iterations for IS:  5000 
+#>     Simulations:
+#>         nb of simulated datasets used for npde:  1000 
+#>         nb of simulated datasets used for VPC:  100 
+#>     Input/output
+#>         save the results to a file:  FALSE 
+#>         save the graphs to files:  FALSE 
+#> ----------------------------------------------------
+#> ----                  Results                   ----
+#> ----------------------------------------------------
+#> -----------------  Fixed effects  ------------------
+#> ----------------------------------------------------
+#>   Parameter    Estimate SE   CV(%)
+#>   parent_0      97.44   <NA> <NA> 
+#>   log_k_parent -56.97   <NA> <NA> 
+#> a a.1           -7.12   <NA> <NA> 
+#> b b.1            0.15   <NA> <NA> 
+#> ----------------------------------------------------
+#> -----------  Variance of random effects  -----------
+#> ----------------------------------------------------
+#>              Parameter           Estimate SE CV(%)
+#> parent_0     omega2.parent_0       2.8    NA NA   
+#> log_k_parent omega2.log_k_parent 750.3    NA NA   
+#> ----------------------------------------------------
+#> ------  Correlation matrix of random effects  ------
+#> ----------------------------------------------------
+#>                     omega2.parent_0 omega2.log_k_parent
+#> omega2.parent_0     1               0                  
+#> omega2.log_k_parent 0               1                  
+#> ----------------------------------------------------
+#> ---------------  Statistical criteria  -------------
+#> ----------------------------------------------------
+#> Likelihood computed by linearisation
+#>       -2LL= 623.7744 
+#>       AIC = 635.7744 
+#>       BIC = 633.431 
+#> 
+#> Likelihood computed by importance sampling
+#>       -2LL= 621.1909 
+#>       AIC = 633.1909 
+#>       BIC = 630.8475 
+#> ----------------------------------------------------
    plot(f_saemix_tc, plot.type = "convergence")
+
    #> Plotting convergence plots
    
 
    
   
    
   
 
 
    
-  
    Site built with pkgdown 1.5.1.
    
+  
    Site built with pkgdown 1.6.1.
    
 
    
 
       
diff --git a/docs/dev/sitemap.xml b/docs/dev/sitemap.xml
index 81368436..e284abf6 100644
--- a/docs/dev/sitemap.xml
+++ b/docs/dev/sitemap.xml
@@ -171,6 +171,9 @@
   
     https://pkgdown.jrwb.de/mkin/reference/residuals.mkinfit.html
   
+  
+    https://pkgdown.jrwb.de/mkin/reference/saemix.html
+  
   
     https://pkgdown.jrwb.de/mkin/reference/schaefer07_complex_case.html
   
diff --git a/man/saemix.Rd b/man/saemix.Rd
index eedbf537..5f96e84f 100644
--- a/man/saemix.Rd
+++ b/man/saemix.Rd
@@ -43,7 +43,7 @@ sfo_sfo <- mkinmod(parent = mkinsub("SFO", "A1"),
 \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,
@@ -51,4 +51,28 @@ saemix_options <- list(seed = 123456,
 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")
+
 }
-- 
cgit v1.2.1