From 6178249bbb5e9de7cb7f34287ee7de28a68fed6c Mon Sep 17 00:00:00 2001 From: Johannes Ranke Date: Wed, 10 Aug 2022 15:38:17 +0200 Subject: Change dev branch used for docs, update static docs --- docs/dev/reference/CAKE_export.html | 282 +++---- docs/dev/reference/D24_2014.html | 313 ++++---- docs/dev/reference/Rplot002.png | Bin 59355 -> 17010 bytes docs/dev/reference/Rplot003.png | Bin 59128 -> 49705 bytes docs/dev/reference/Rplot004.png | Bin 49360 -> 58906 bytes docs/dev/reference/Rplot005.png | Bin 59600 -> 22787 bytes docs/dev/reference/Rplot006.png | Bin 22144 -> 24295 bytes docs/dev/reference/Rplot007.png | Bin 24992 -> 25204 bytes docs/dev/reference/convergence.html | 163 ++++ docs/dev/reference/create_deg_func.html | 228 +++--- docs/dev/reference/dimethenamid_2018-1.png | Bin 263830 -> 248143 bytes docs/dev/reference/dimethenamid_2018-2.png | Bin 246327 -> 246117 bytes docs/dev/reference/dimethenamid_2018-3.png | Bin 0 -> 247960 bytes docs/dev/reference/dimethenamid_2018.html | 445 ++++------- docs/dev/reference/endpoints.html | 40 +- docs/dev/reference/illparms.html | 207 +++++ docs/dev/reference/index.html | 61 +- docs/dev/reference/intervals.saem.mmkin.html | 24 +- docs/dev/reference/mean_degparms.html | 34 +- docs/dev/reference/mhmkin.html | 213 ++++++ docs/dev/reference/mkinds.html | 277 +++---- docs/dev/reference/mkindsg.html | 689 ++++++++--------- docs/dev/reference/mkinfit-1.png | Bin 65835 -> 66250 bytes docs/dev/reference/mkinfit.html | 1050 +++++++++++++------------- docs/dev/reference/nlme.mmkin-1.png | Bin 124929 -> 124053 bytes docs/dev/reference/nlme.mmkin-2.png | Bin 169786 -> 168550 bytes docs/dev/reference/nlme.mmkin-3.png | Bin 173135 -> 172051 bytes docs/dev/reference/nlme.mmkin.html | 239 +++--- docs/dev/reference/plot.mixed.mmkin-2.png | Bin 174398 -> 173322 bytes docs/dev/reference/plot.mixed.mmkin-3.png | Bin 173794 -> 172513 bytes docs/dev/reference/plot.mixed.mmkin-4.png | Bin 176972 -> 175535 bytes docs/dev/reference/plot.mixed.mmkin.html | 171 +++-- docs/dev/reference/plot.mkinfit-1.png | Bin 53151 -> 53304 bytes docs/dev/reference/plot.mkinfit-2.png | Bin 73254 -> 73194 bytes docs/dev/reference/plot.mkinfit-3.png | Bin 67810 -> 68138 bytes docs/dev/reference/plot.mkinfit-4.png | Bin 72295 -> 72694 bytes docs/dev/reference/plot.mkinfit-5.png | Bin 66560 -> 67141 bytes docs/dev/reference/plot.mkinfit-6.png | Bin 72717 -> 72904 bytes docs/dev/reference/plot.mkinfit-7.png | Bin 73553 -> 74323 bytes docs/dev/reference/plot.mkinfit.html | 532 ++++++------- docs/dev/reference/reexports.html | 19 +- docs/dev/reference/saem-1.png | Bin 47337 -> 46419 bytes docs/dev/reference/saem-2.png | Bin 48793 -> 49282 bytes docs/dev/reference/saem-3.png | Bin 82192 -> 128227 bytes docs/dev/reference/saem-4.png | Bin 128209 -> 171244 bytes docs/dev/reference/saem.html | 613 +++++++++++---- docs/dev/reference/summary.mkinfit.html | 439 +++++------ docs/dev/reference/summary.mmkin.html | 173 +++++ docs/dev/reference/summary.nlme.mmkin.html | 658 ++++++++-------- docs/dev/reference/summary.saem.mmkin.html | 454 +++++++++-- 50 files changed, 4194 insertions(+), 3130 deletions(-) create mode 100644 docs/dev/reference/convergence.html create mode 100644 docs/dev/reference/dimethenamid_2018-3.png create mode 100644 docs/dev/reference/illparms.html create mode 100644 docs/dev/reference/mhmkin.html create mode 100644 docs/dev/reference/summary.mmkin.html (limited to 'docs/dev/reference') diff --git a/docs/dev/reference/CAKE_export.html b/docs/dev/reference/CAKE_export.html index e187772f..47efd056 100644 --- a/docs/dev/reference/CAKE_export.html +++ b/docs/dev/reference/CAKE_export.html @@ -1,68 +1,13 @@ - - - - - - - -Export a list of datasets format to a CAKE study file — CAKE_export • mkin - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -Export a list of datasets format to a CAKE study file — CAKE_export • mkin - + + - - - -
-
- -
- -
+

Export a list of datasets format to a CAKE study file

- Source: R/CAKE_export.R + Source: R/CAKE_export.R
@@ -149,116 +88,115 @@ specified as well." /> specified as well.

- 
CAKE_export(
-  ds,
-  map = c(parent = "Parent"),
-  links = NA,
-  filename = "CAKE_export.csf",
-  path = ".",
-  overwrite = FALSE,
-  study = "Codlemone aerobic soil degradation",
-  description = "",
-  time_unit = "days",
-  res_unit = "% AR",
-  comment = "Created using mkin::CAKE_export",
-  date = Sys.Date(),
-  optimiser = "IRLS"
-)
- -

Arguments

- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
ds

A named list of datasets in long format as compatible with -mkinfit.

map

A character vector with CAKE compartment names (Parent, A1, ...), -named with the names used in the list of datasets.

links

An optional character vector of target compartments, named with +

+ 
CAKE_export(
+  ds,
+  map = c(parent = "Parent"),
+  links = NA,
+  filename = "CAKE_export.csf",
+  path = ".",
+  overwrite = FALSE,
+  study = "Degradinol aerobic soil degradation",
+  description = "",
+  time_unit = "days",
+  res_unit = "% AR",
+  comment = "",
+  date = Sys.Date(),
+  optimiser = "IRLS"
+)
+
+ +
+

Arguments

+
ds
+

A named list of datasets in long format as compatible with +mkinfit.

+ + +
map
+

A character vector with CAKE compartment names (Parent, A1, ...), +named with the names used in the list of datasets.

+ + +
links
+

An optional character vector of target compartments, named with the names of the source compartments. In order to make this easier, the -names are used as in the datasets supplied.

filename

Where to write the result. Should end in .csf in order to be -compatible with CAKE.

path

An optional path to the output file.

overwrite

If TRUE, existing files are overwritten.

study

The name of the study.

description

An optional description.

time_unit

The time unit for the residue data.

res_unit

The unit used for the residues.

comment

An optional comment.

date

The date of file creation.

optimiser

Can be OLS or IRLS.

- -

Value

- -

The function is called for its side effect.

-

Author

+names are used as in the datasets supplied.

+ + +
filename
+

Where to write the result. Should end in .csf in order to be +compatible with CAKE.

+ + +
path
+

An optional path to the output file.

+ + +
overwrite
+

If TRUE, existing files are overwritten.

+ + +
study
+

The name of the study.

+ + +
description
+

An optional description.

+ +
time_unit
+

The time unit for the residue data.

+ + +
res_unit
+

The unit used for the residues.

+ + +
comment
+

An optional comment.

+ + +
date
+

The date of file creation.

+ + +
optimiser
+

Can be OLS or IRLS.

+ +
+
+

Value

+ + +

The function is called for its side effect.

+
+
+

Author

Johannes Ranke

+
+
- - - + + diff --git a/docs/dev/reference/D24_2014.html b/docs/dev/reference/D24_2014.html index 9ecb6330..5cf7604c 100644 --- a/docs/dev/reference/D24_2014.html +++ b/docs/dev/reference/D24_2014.html @@ -1,72 +1,17 @@ - - - - - - - -Aerobic soil degradation data on 2,4-D from the EU assessment in 2014 — D24_2014 • mkin - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -Aerobic soil degradation data on 2,4-D from the EU assessment in 2014 — D24_2014 • mkin - - - - - - - - - - + + - - - -
-
- -
- -
+

Aerobic soil degradation data on 2,4-D from the EU assessment in 2014

- Source: R/D24_2014.R + Source: R/D24_2014.R
@@ -157,20 +96,23 @@ context of pesticide registrations, as the use of the data may be constrained by data protection regulations.

- 
D24_2014
- - -

Format

- -

An mkindsg object grouping five datasets

-

Source

+
+ 
D24_2014
+
+
+

Format

+

An mkindsg object grouping five datasets

+
+
+

Source

Hellenic Ministry of Rural Development and Agriculture (2014) Final addendum to the Renewal Assessment Report - public version - 2,4-D Volume 3 Annex B.8 Fate and behaviour in the environment -http://registerofquestions.efsa.europa.eu/roqFrontend/outputLoader?output=ON-3812

-

Details

- +https://open.efsa.europa.eu/study-inventory/EFSA-Q-2013-00811

+
+
+

Details

Data for the first dataset are from p. 685. Data for the other four datasets were used in the preprocessed versions given in the kinetics section (p. 761ff.), with the exception of residues smaller than 1 for DCP @@ -178,114 +120,119 @@ in the soil from Site I2, where the values given on p. 694 were used.

The R code used to create this data object is installed with this package in the 'dataset_generation' directory. In the code, page numbers are given for specific pieces of information in the comments.

+
-

Examples

- 
print(D24_2014)
-
#> <mkindsg> holding 5 mkinds objects
-#> Title $title:  Aerobic soil degradation data on 2,4-D from the EU assessment in 2014 
-#> Occurrence of observed compounds $observed_n:
-#> D24 DCP DCA 
-#>   5   4   4 
-#> Time normalisation factors $f_time_norm:
-#> [1] 1.6062378 0.7118732 0.7156063 0.7156063 0.8977124
-#> Meta information $meta:
-#>                                 study usda_soil_type study_moisture_ref_type
-#> Mississippi                Cohen 1991      Silt loam                    <NA>
-#> Fayette     Liu and Adelfinskaya 2011      Silt loam                     pF1
-#> RefSol 03-G Liu and Adelfinskaya 2011           Loam                     pF1
-#> Site E1     Liu and Adelfinskaya 2011           Loam                     pF1
-#> Site I2     Liu and Adelfinskaya 2011     Loamy sand                     pF1
-#>             rel_moisture temperature
-#> Mississippi           NA          25
-#> Fayette              0.5          20
-#> RefSol 03-G          0.5          20
-#> Site E1              0.5          20
-#> Site I2              0.5          20
# \dontrun{
-print(D24_2014$ds[[1]], data = TRUE)
-
#> <mkinds> with $title:  Mississippi 
-#> Observed compounds $observed:  D24 
-#> Sampling times $sampling_times:
-#> 0, 2, 4, 7, 15, 24, 35, 56, 71, 114, 183, 273, 365 
-#> With a maximum of  1  replicates
-#>    time  D24
-#> 1     0 96.8
-#> 2     2 81.0
-#> 3     4 81.7
-#> 4     7 88.2
-#> 5    15 66.3
-#> 6    24 72.9
-#> 7    35 62.6
-#> 8    56 54.6
-#> 9    71 35.2
-#> 10  114 18.0
-#> 11  183 11.3
-#> 12  273  9.9
-#> 13  365  6.3
m_D24 = mkinmod(D24 = mkinsub("SFO", to = "DCP"),
-  DCP = mkinsub("SFO", to = "DCA"),
-  DCA = mkinsub("SFO"))
-
#> Temporary DLL for differentials generated and loaded
print(m_D24)
-
#> <mkinmod> model generated with
-#> Use of formation fractions $use_of_ff: max 
-#> Specification $spec:
-#> $D24
-#> $type: SFO; $to: DCP; $sink: TRUE
-#> $DCP
-#> $type: SFO; $to: DCA; $sink: TRUE
-#> $DCA
-#> $type: SFO; $sink: TRUE
-#> Coefficient matrix $coefmat available
-#> Compiled model $cf available
-#> Differential equations:
-#> d_D24/dt = - k_D24 * D24
-#> d_DCP/dt = + f_D24_to_DCP * k_D24 * D24 - k_DCP * DCP
-#> d_DCA/dt = + f_DCP_to_DCA * k_DCP * DCP - k_DCA * DCA
m_D24_2 = mkinmod(D24 = mkinsub("DFOP", to = "DCP"),
-  DCP = mkinsub("SFO", to = "DCA"),
-  DCA = mkinsub("SFO"))
-
#> Temporary DLL for differentials generated and loaded
print(m_D24_2)
-
#> <mkinmod> model generated with
-#> Use of formation fractions $use_of_ff: max 
-#> Specification $spec:
-#> $D24
-#> $type: DFOP; $to: DCP; $sink: TRUE
-#> $DCP
-#> $type: SFO; $to: DCA; $sink: TRUE
-#> $DCA
-#> $type: SFO; $sink: TRUE
-#> Compiled model $cf available
-#> Differential equations:
-#> d_D24/dt = - ((k1 * g * exp(-k1 * time) + k2 * (1 - g) * exp(-k2 *
-#>            time)) / (g * exp(-k1 * time) + (1 - g) * exp(-k2 * time)))
-#>            * D24
-#> d_DCP/dt = + f_D24_to_DCP * ((k1 * g * exp(-k1 * time) + k2 * (1 - g) *
-#>            exp(-k2 * time)) / (g * exp(-k1 * time) + (1 - g) * exp(-k2
-#>            * time))) * D24 - k_DCP * DCP
-#> d_DCA/dt = + f_DCP_to_DCA * k_DCP * DCP - k_DCA * DCA
# }
-
+
+

Examples

+ 
print(D24_2014)
+#> <mkindsg> holding 5 mkinds objects
+#> Title $title:  Aerobic soil degradation data on 2,4-D from the EU assessment in 2014 
+#> Occurrence of observed compounds $observed_n:
+#> D24 DCP DCA 
+#>   5   4   4 
+#> Time normalisation factors $f_time_norm:
+#> [1] 1.6062378 0.7118732 0.7156063 0.7156063 0.8977124
+#> Meta information $meta:
+#>                                 study usda_soil_type study_moisture_ref_type
+#> Mississippi                Cohen 1991      Silt loam                    <NA>
+#> Fayette     Liu and Adelfinskaya 2011      Silt loam                     pF1
+#> RefSol 03-G Liu and Adelfinskaya 2011           Loam                     pF1
+#> Site E1     Liu and Adelfinskaya 2011           Loam                     pF1
+#> Site I2     Liu and Adelfinskaya 2011     Loamy sand                     pF1
+#>             rel_moisture temperature
+#> Mississippi           NA          25
+#> Fayette              0.5          20
+#> RefSol 03-G          0.5          20
+#> Site E1              0.5          20
+#> Site I2              0.5          20
+# \dontrun{
+print(D24_2014$ds[[1]], data = TRUE)
+#> <mkinds> with $title:  Mississippi 
+#> Observed compounds $observed:  D24 
+#> Sampling times $sampling_times:
+#> 0, 2, 4, 7, 15, 24, 35, 56, 71, 114, 183, 273, 365 
+#> With a maximum of  1  replicates
+#>    time  D24
+#> 1     0 96.8
+#> 2     2 81.0
+#> 3     4 81.7
+#> 4     7 88.2
+#> 5    15 66.3
+#> 6    24 72.9
+#> 7    35 62.6
+#> 8    56 54.6
+#> 9    71 35.2
+#> 10  114 18.0
+#> 11  183 11.3
+#> 12  273  9.9
+#> 13  365  6.3
+m_D24 = mkinmod(D24 = mkinsub("SFO", to = "DCP"),
+  DCP = mkinsub("SFO", to = "DCA"),
+  DCA = mkinsub("SFO"))
+#> Temporary DLL for differentials generated and loaded
+print(m_D24)
+#> <mkinmod> model generated with
+#> Use of formation fractions $use_of_ff: max 
+#> Specification $spec:
+#> $D24
+#> $type: SFO; $to: DCP; $sink: TRUE
+#> $DCP
+#> $type: SFO; $to: DCA; $sink: TRUE
+#> $DCA
+#> $type: SFO; $sink: TRUE
+#> Coefficient matrix $coefmat available
+#> Compiled model $cf available
+#> Differential equations:
+#> d_D24/dt = - k_D24 * D24
+#> d_DCP/dt = + f_D24_to_DCP * k_D24 * D24 - k_DCP * DCP
+#> d_DCA/dt = + f_DCP_to_DCA * k_DCP * DCP - k_DCA * DCA
+m_D24_2 = mkinmod(D24 = mkinsub("DFOP", to = "DCP"),
+  DCP = mkinsub("SFO", to = "DCA"),
+  DCA = mkinsub("SFO"))
+#> Temporary DLL for differentials generated and loaded
+print(m_D24_2)
+#> <mkinmod> model generated with
+#> Use of formation fractions $use_of_ff: max 
+#> Specification $spec:
+#> $D24
+#> $type: DFOP; $to: DCP; $sink: TRUE
+#> $DCP
+#> $type: SFO; $to: DCA; $sink: TRUE
+#> $DCA
+#> $type: SFO; $sink: TRUE
+#> Compiled model $cf available
+#> Differential equations:
+#> d_D24/dt = - ((k1 * g * exp(-k1 * time) + k2 * (1 - g) * exp(-k2 *
+#>            time)) / (g * exp(-k1 * time) + (1 - g) * exp(-k2 * time)))
+#>            * D24
+#> d_DCP/dt = + f_D24_to_DCP * ((k1 * g * exp(-k1 * time) + k2 * (1 - g) *
+#>            exp(-k2 * time)) / (g * exp(-k1 * time) + (1 - g) * exp(-k2
+#>            * time))) * D24 - k_DCP * DCP
+#> d_DCA/dt = + f_DCP_to_DCA * k_DCP * DCP - k_DCA * DCA
+# }
+
+
+
- - - + + diff --git a/docs/dev/reference/Rplot002.png b/docs/dev/reference/Rplot002.png index e40f40a0..f06a860e 100644 Binary files a/docs/dev/reference/Rplot002.png and b/docs/dev/reference/Rplot002.png differ diff --git a/docs/dev/reference/Rplot003.png b/docs/dev/reference/Rplot003.png index fa29fc43..1af5d4b4 100644 Binary files a/docs/dev/reference/Rplot003.png and b/docs/dev/reference/Rplot003.png differ diff --git a/docs/dev/reference/Rplot004.png b/docs/dev/reference/Rplot004.png index c279f831..12d337a4 100644 Binary files a/docs/dev/reference/Rplot004.png and b/docs/dev/reference/Rplot004.png differ diff --git a/docs/dev/reference/Rplot005.png b/docs/dev/reference/Rplot005.png index 92c7cc2d..cb419daa 100644 Binary files a/docs/dev/reference/Rplot005.png and b/docs/dev/reference/Rplot005.png differ diff --git a/docs/dev/reference/Rplot006.png b/docs/dev/reference/Rplot006.png index 4c728f4e..48f5bbd8 100644 Binary files a/docs/dev/reference/Rplot006.png and b/docs/dev/reference/Rplot006.png differ diff --git a/docs/dev/reference/Rplot007.png b/docs/dev/reference/Rplot007.png index 10b7455a..21a6ea76 100644 Binary files a/docs/dev/reference/Rplot007.png and b/docs/dev/reference/Rplot007.png differ diff --git a/docs/dev/reference/convergence.html b/docs/dev/reference/convergence.html new file mode 100644 index 00000000..e9fac3bb --- /dev/null +++ b/docs/dev/reference/convergence.html @@ -0,0 +1,163 @@ + +Method to get convergence information — convergence • mkin + + +
+
+ + + +
+
+
+

Method to get convergence information

+ Source: R/convergence.R + +
+ +
+

Method to get convergence information

+
+ +
+ 
convergence(object, ...)
+
+# S3 method for mmkin
+convergence(object, ...)
+
+# S3 method for convergence.mmkin
+print(x, ...)
+
+ +
+

Arguments

+
object
+

The object to investigate

+ + +
...
+

For potential future extensions

+ + +
x
+

The object to be printed

+ +
+
+

Value

+ + +

For mkinfit objects, a character vector containing +For mmkin objects, an object of class 'convergence.mmkin' with a +suitable printing method.

+
+ +
+

Examples

+ 
# \dontrun{
+fits <- mmkin(
+  c("SFO", "FOMC"),
+  list("FOCUS A" = FOCUS_2006_A,
+       "FOCUS B" = FOCUS_2006_C),
+  quiet = TRUE)
+convergence(fits)
+#>       dataset
+#> model  FOCUS A FOCUS B
+#>   SFO  OK      OK     
+#>   FOMC OK      OK     
+#> 
+#> OK: No warnings
+# }
+
+
+
+ +
+ + +
+

Developed by Johannes Ranke.

+
+ +
+

Site built with pkgdown 2.0.6.

+
+ +
+ + + + + + + + diff --git a/docs/dev/reference/create_deg_func.html b/docs/dev/reference/create_deg_func.html index 65a682bb..5d5870fe 100644 --- a/docs/dev/reference/create_deg_func.html +++ b/docs/dev/reference/create_deg_func.html @@ -1,67 +1,12 @@ - - - - - - - -Create degradation functions for known analytical solutions — create_deg_func • mkin - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -Create degradation functions for known analytical solutions — create_deg_func • mkin - + + - - - -
-
- -
- -
+

Create degradation functions for known analytical solutions

- Source: R/create_deg_func.R + Source: R/create_deg_func.R
@@ -147,83 +86,88 @@

Create degradation functions for known analytical solutions

- 
create_deg_func(spec, use_of_ff = c("min", "max"))
- -

Arguments

- - - - - - - - - - -
spec

List of model specifications as contained in mkinmod objects

use_of_ff

Minimum or maximum use of formation fractions

- -

Value

- -

Degradation function to be attached to mkinmod objects

- -

Examples

- 

-SFO_SFO <- mkinmod(
-  parent = mkinsub("SFO", "m1"),
-  m1 = mkinsub("SFO"))
-
#> Temporary DLL for differentials generated and loaded
FOCUS_D <- subset(FOCUS_2006_D, value != 0) # to avoid warnings
-fit_1 <- mkinfit(SFO_SFO, FOCUS_D, solution_type = "analytical", quiet = TRUE)
-fit_2 <- mkinfit(SFO_SFO, FOCUS_D, solution_type = "deSolve", quiet = TRUE)
-# \dontrun{
-if (require(rbenchmark))
-  benchmark(
-    analytical = mkinfit(SFO_SFO, FOCUS_D, solution_type = "analytical", quiet = TRUE),
-    deSolve = mkinfit(SFO_SFO, FOCUS_D, solution_type = "deSolve", quiet = TRUE),
-    replications = 2)
-
#> Loading required package: rbenchmark
#>         test replications elapsed relative user.self sys.self user.child
-#> 1 analytical            2   0.396     1.00     0.396        0          0
-#> 2    deSolve            2   0.709     1.79     0.707        0          0
-#>   sys.child
-#> 1         0
-#> 2         0
  DFOP_SFO <- mkinmod(
-    parent = mkinsub("DFOP", "m1"),
-    m1 = mkinsub("SFO"))
-
#> Temporary DLL for differentials generated and loaded
  benchmark(
-    analytical = mkinfit(DFOP_SFO, FOCUS_D, solution_type = "analytical", quiet = TRUE),
-    deSolve = mkinfit(DFOP_SFO, FOCUS_D, solution_type = "deSolve", quiet = TRUE),
-    replications = 2)
-
#>         test replications elapsed relative user.self sys.self user.child
-#> 1 analytical            2   0.844    1.000     0.844        0          0
-#> 2    deSolve            2   1.624    1.924     1.624        0          0
-#>   sys.child
-#> 1         0
-#> 2         0
# }
-
+
+ 
create_deg_func(spec, use_of_ff = c("min", "max"))
+
+ +
+

Arguments

+
spec
+

List of model specifications as contained in mkinmod objects

+ + +
use_of_ff
+

Minimum or maximum use of formation fractions

+ +
+
+

Value

+ + +

Degradation function to be attached to mkinmod objects

+
+ +
+

Examples

+ 

+SFO_SFO <- mkinmod(
+  parent = mkinsub("SFO", "m1"),
+  m1 = mkinsub("SFO"))
+#> Temporary DLL for differentials generated and loaded
+FOCUS_D <- subset(FOCUS_2006_D, value != 0) # to avoid warnings
+fit_1 <- mkinfit(SFO_SFO, FOCUS_D, solution_type = "analytical", quiet = TRUE)
+# \dontrun{
+fit_2 <- mkinfit(SFO_SFO, FOCUS_D, solution_type = "deSolve", quiet = TRUE)
+if (require(rbenchmark))
+  benchmark(
+    analytical = mkinfit(SFO_SFO, FOCUS_D, solution_type = "analytical", quiet = TRUE),
+    deSolve = mkinfit(SFO_SFO, FOCUS_D, solution_type = "deSolve", quiet = TRUE),
+    replications = 2)
+#> Loading required package: rbenchmark
+#>         test replications elapsed relative user.self sys.self user.child
+#> 1 analytical            2   0.401    1.000     0.401        0          0
+#> 2    deSolve            2   0.654    1.631     0.654        0          0
+#>   sys.child
+#> 1         0
+#> 2         0
+  DFOP_SFO <- mkinmod(
+    parent = mkinsub("DFOP", "m1"),
+    m1 = mkinsub("SFO"))
+#> Temporary DLL for differentials generated and loaded
+  benchmark(
+    analytical = mkinfit(DFOP_SFO, FOCUS_D, solution_type = "analytical", quiet = TRUE),
+    deSolve = mkinfit(DFOP_SFO, FOCUS_D, solution_type = "deSolve", quiet = TRUE),
+    replications = 2)
+#>         test replications elapsed relative user.self sys.self user.child
+#> 1 analytical            2   0.843    1.000     0.843        0          0
+#> 2    deSolve            2   1.445    1.714     1.445        0          0
+#>   sys.child
+#> 1         0
+#> 2         0
+# }
+
+
+
- - - + + diff --git a/docs/dev/reference/dimethenamid_2018-1.png b/docs/dev/reference/dimethenamid_2018-1.png index b8c5355f..4300b0c0 100644 Binary files a/docs/dev/reference/dimethenamid_2018-1.png and b/docs/dev/reference/dimethenamid_2018-1.png differ diff --git a/docs/dev/reference/dimethenamid_2018-2.png b/docs/dev/reference/dimethenamid_2018-2.png index 3b8a123b..36db063c 100644 Binary files a/docs/dev/reference/dimethenamid_2018-2.png and b/docs/dev/reference/dimethenamid_2018-2.png differ diff --git a/docs/dev/reference/dimethenamid_2018-3.png b/docs/dev/reference/dimethenamid_2018-3.png new file mode 100644 index 00000000..7c876208 Binary files /dev/null and b/docs/dev/reference/dimethenamid_2018-3.png differ diff --git a/docs/dev/reference/dimethenamid_2018.html b/docs/dev/reference/dimethenamid_2018.html index 5fb94988..7a356284 100644 --- a/docs/dev/reference/dimethenamid_2018.html +++ b/docs/dev/reference/dimethenamid_2018.html @@ -22,7 +22,7 @@ constrained by data protection regulations."> mkin - 1.1.0 + 1.1.2 @@ -31,7 +31,7 @@ constrained by data protection regulations.">Functions and data
  • -
    • +
  • + Example evaluations of dimethenamid data from 2018 with nonlinear mixed-effects models +
  • Example evaluation of FOCUS Example Dataset Z
    • @@ -94,7 +97,7 @@ constrained by data protection regulations.

    - 
    dimethenamid_2018
    + 
    dimethenamid_2018
    @@ -117,7 +120,7 @@ specific pieces of information in the comments.

    Examples

    - 
    print(dimethenamid_2018)
+    
    print(dimethenamid_2018)
 #> <mkindsg> holding 7 mkinds objects
 #> Title $title:  Aerobic soil degradation data on dimethenamid-P from the EU assessment in 2018 
 #> Occurrence of observed compounds $observed_n:
@@ -142,296 +145,168 @@ specific pieces of information in the comments.
    
 #> Flaach                   NA          20
 #> BBA 2.2                  NA          20
 #> BBA 2.3                  NA          20
-dmta_ds <- lapply(1:7, function(i) {
-  ds_i <- dimethenamid_2018$ds[[i]]$data
-  ds_i[ds_i$name == "DMTAP", "name"] <-  "DMTA"
-  ds_i$time <- ds_i$time * dimethenamid_2018$f_time_norm[i]
-  ds_i
-})
-names(dmta_ds) <- sapply(dimethenamid_2018$ds, function(ds) ds$title)
-dmta_ds[["Elliot"]] <- rbind(dmta_ds[["Elliot 1"]], dmta_ds[["Elliot 2"]])
-dmta_ds[["Elliot 1"]] <- NULL
-dmta_ds[["Elliot 2"]] <- NULL
-# \dontrun{
-dfop_sfo3_plus <- mkinmod(
-  DMTA = mkinsub("DFOP", c("M23", "M27", "M31")),
-  M23 = mkinsub("SFO"),
-  M27 = mkinsub("SFO"),
-  M31 = mkinsub("SFO", "M27", sink = FALSE),
-  quiet = TRUE
-)
-f_dmta_mkin_tc <- mmkin(
-  list("DFOP-SFO3+" = dfop_sfo3_plus),
-  dmta_ds, quiet = TRUE, error_model = "tc")
-nlmixr_model(f_dmta_mkin_tc)
-#> With est = 'saem', a different error model is required for each observed variableChanging the error model to 'obs_tc' (Two-component error for each observed variable)
-#> Warning: number of items to replace is not a multiple of replacement length
-#> function () 
-#> {
-#>     ini({
-#>         DMTA_0 = 99
-#>         eta.DMTA_0 ~ 2.3
-#>         log_k_M23 = -3.9
-#>         eta.log_k_M23 ~ 0.55
-#>         log_k_M27 = -4.3
-#>         eta.log_k_M27 ~ 0.86
-#>         log_k_M31 = -4.2
-#>         eta.log_k_M31 ~ 0.75
-#>         log_k1 = -2.2
-#>         eta.log_k1 ~ 0.9
-#>         log_k2 = -3.8
-#>         eta.log_k2 ~ 1.6
-#>         g_qlogis = 0.44
-#>         eta.g_qlogis ~ 3.1
-#>         f_DMTA_tffm0_1_qlogis = -2.1
-#>         eta.f_DMTA_tffm0_1_qlogis ~ 0.3
-#>         f_DMTA_tffm0_2_qlogis = -2.2
-#>         eta.f_DMTA_tffm0_2_qlogis ~ 0.3
-#>         f_DMTA_tffm0_3_qlogis = -2.1
-#>         eta.f_DMTA_tffm0_3_qlogis ~ 0.3
-#>         sigma_low_DMTA = 0.7
-#>         rsd_high_DMTA = 0.026
-#>         sigma_low_M23 = 0.7
-#>         rsd_high_M23 = 0.026
-#>         sigma_low_M27 = 0.7
-#>         rsd_high_M27 = 0.026
-#>         sigma_low_M31 = 0.7
-#>         rsd_high_M31 = 0.026
-#>     })
-#>     model({
-#>         DMTA_0_model = DMTA_0 + eta.DMTA_0
-#>         DMTA(0) = DMTA_0_model
-#>         k_M23 = exp(log_k_M23 + eta.log_k_M23)
-#>         k_M27 = exp(log_k_M27 + eta.log_k_M27)
-#>         k_M31 = exp(log_k_M31 + eta.log_k_M31)
-#>         k1 = exp(log_k1 + eta.log_k1)
-#>         k2 = exp(log_k2 + eta.log_k2)
-#>         g = expit(g_qlogis + eta.g_qlogis)
-#>         f_DMTA_to_M23 = expit(f_DMTA_tffm0_1_qlogis + eta.f_DMTA_tffm0_1_qlogis)
-#>         f_DMTA_to_M23 = expit(f_DMTA_tffm0_2_qlogis + eta.f_DMTA_tffm0_2_qlogis)
-#>         f_DMTA_to_M23 = expit(f_DMTA_tffm0_3_qlogis + eta.f_DMTA_tffm0_3_qlogis)
-#>         f_DMTA_to_M23 = f_DMTA_tffm0_1
-#>         f_DMTA_to_M27 = f_DMTA_tffm0_2 * (1 - f_DMTA_tffm0_1)
-#>         f_DMTA_to_M31 = f_DMTA_tffm0_3 * (1 - f_DMTA_tffm0_2) * 
-#>             (1 - f_DMTA_tffm0_1)
-#>         d/dt(DMTA) = -((k1 * g * exp(-k1 * time) + k2 * (1 - 
-#>             g) * exp(-k2 * time))/(g * exp(-k1 * time) + (1 - 
-#>             g) * exp(-k2 * time))) * DMTA
-#>         d/dt(M23) = +f_DMTA_to_M23 * ((k1 * g * exp(-k1 * time) + 
-#>             k2 * (1 - g) * exp(-k2 * time))/(g * exp(-k1 * time) + 
-#>             (1 - g) * exp(-k2 * time))) * DMTA - k_M23 * M23
-#>         d/dt(M27) = +f_DMTA_to_M27 * ((k1 * g * exp(-k1 * time) + 
-#>             k2 * (1 - g) * exp(-k2 * time))/(g * exp(-k1 * time) + 
-#>             (1 - g) * exp(-k2 * time))) * DMTA - k_M27 * M27 + 
-#>             k_M31 * M31
-#>         d/dt(M31) = +f_DMTA_to_M31 * ((k1 * g * exp(-k1 * time) + 
-#>             k2 * (1 - g) * exp(-k2 * time))/(g * exp(-k1 * time) + 
-#>             (1 - g) * exp(-k2 * time))) * DMTA - k_M31 * M31
-#>         DMTA ~ add(sigma_low_DMTA) + prop(rsd_high_DMTA)
-#>         M23 ~ add(sigma_low_M23) + prop(rsd_high_M23)
-#>         M27 ~ add(sigma_low_M27) + prop(rsd_high_M27)
-#>         M31 ~ add(sigma_low_M31) + prop(rsd_high_M31)
-#>     })
-#> }
-#> <environment: 0x55555fca3790>
-# The focei fit takes about four minutes on my system
-system.time(
-  f_dmta_nlmixr_focei <- nlmixr(f_dmta_mkin_tc, est = "focei",
-    control = nlmixr::foceiControl(print = 500))
-)
-#> Warning: number of items to replace is not a multiple of replacement length
-#> ℹ parameter labels from comments are typically ignored in non-interactive mode
-#> ℹ Need to run with the source intact to parse comments
-#> → creating full model...
-#> → pruning branches (`if`/`else`)...
-#> ✔ done
-#> → loading into symengine environment...
-#> ✔ done
-#> → creating full model...
-#> → pruning branches (`if`/`else`)...
-#> ✔ done
-#> → loading into symengine environment...
-#> ✔ done
-#> → calculate jacobian
-#> [====|====|====|====|====|====|====|====|====|====] 0:00:01 
+dmta_ds <- lapply(1:7, function(i) {
+  ds_i <- dimethenamid_2018$ds[[i]]$data
+  ds_i[ds_i$name == "DMTAP", "name"] <-  "DMTA"
+  ds_i$time <- ds_i$time * dimethenamid_2018$f_time_norm[i]
+  ds_i
+})
+names(dmta_ds) <- sapply(dimethenamid_2018$ds, function(ds) ds$title)
+dmta_ds[["Elliot"]] <- rbind(dmta_ds[["Elliot 1"]], dmta_ds[["Elliot 2"]])
+dmta_ds[["Elliot 1"]] <- NULL
+dmta_ds[["Elliot 2"]] <- NULL
+# \dontrun{
+# We don't use DFOP for the parent compound, as this gives numerical
+# instabilities in the fits
+sfo_sfo3p <- mkinmod(
+ DMTA = mkinsub("SFO", c("M23", "M27", "M31")),
+ M23 = mkinsub("SFO"),
+ M27 = mkinsub("SFO"),
+ M31 = mkinsub("SFO", "M27", sink = FALSE),
+ quiet = TRUE
+)
+dmta_sfo_sfo3p_tc <- mmkin(list("SFO-SFO3+" = sfo_sfo3p),
+  dmta_ds, error_model = "tc", quiet = TRUE)
+print(dmta_sfo_sfo3p_tc)
+#> <mmkin> object
+#> Status of individual fits:
+#> 
+#>            dataset
+#> model       Calke Borstel Flaach BBA 2.2 BBA 2.3 Elliot
+#>   SFO-SFO3+ OK    OK      OK     OK      OK      OK    
+#> 
+#> OK: No warnings
+# The default (test_log_parms = FALSE) gives an undue
+# influence of ill-defined rate constants that have
+# extremely small values:
+plot(mixed(dmta_sfo_sfo3p_tc), test_log_parms = FALSE)
+
+# If we disregards ill-defined rate constants, the results
+# look more plausible, but the truth is likely to be in
+# between these variants
+plot(mixed(dmta_sfo_sfo3p_tc), test_log_parms = TRUE)
+
+# We can also specify a default value for the failing
+# log parameters, to mimic FOCUS guidance
+plot(mixed(dmta_sfo_sfo3p_tc), test_log_parms = TRUE,
+  default_log_parms = log(2)/1000)
+
+# As these attempts are not satisfying, we use nonlinear mixed-effects models
+# f_dmta_nlme_tc <- nlme(dmta_sfo_sfo3p_tc)
+# nlme reaches maxIter = 50 without convergence
+f_dmta_saem_tc <- saem(dmta_sfo_sfo3p_tc)
+# I am commenting out the convergence plot as rendering them
+# with pkgdown fails (at least without further tweaks to the
+# graphics device used)
+#saemix::plot(f_dmta_saem_tc$so, plot.type = "convergence")
+summary(f_dmta_saem_tc)
+#> saemix version used for fitting:      3.1 
+#> mkin version used for pre-fitting:  1.1.2 
+#> R version used for fitting:         4.2.1 
+#> Date of fit:     Wed Aug 10 15:24:12 2022 
+#> Date of summary: Wed Aug 10 15:24:12 2022 
+#> 
+#> Equations:
+#> d_DMTA/dt = - k_DMTA * DMTA
+#> d_M23/dt = + f_DMTA_to_M23 * k_DMTA * DMTA - k_M23 * M23
+#> d_M27/dt = + f_DMTA_to_M27 * k_DMTA * DMTA - k_M27 * M27 + k_M31 * M31
+#> d_M31/dt = + f_DMTA_to_M31 * k_DMTA * DMTA - k_M31 * M31
+#> 
+#> Data:
+#> 563 observations of 4 variable(s) grouped in 6 datasets
+#> 
+#> Model predictions using solution type deSolve 
+#> 
+#> Fitted in 791.863 s
+#> Using 300, 100 iterations and 9 chains
+#> 
+#> Variance model: Two-component variance function 
+#> 
+#> Mean of starting values for individual parameters:
+#>       DMTA_0   log_k_DMTA    log_k_M23    log_k_M27    log_k_M31 f_DMTA_ilr_1 
+#>      95.5662      -2.9048      -3.8130      -4.1600      -4.1486       0.1341 
+#> f_DMTA_ilr_2 f_DMTA_ilr_3 
+#>       0.1385      -1.6700 
 #> 
-#> → calculate sensitivities
-#> [====|====|====|====|====|====|====|====|====|====] 0:00:03 
+#> Fixed degradation parameter values:
+#> None
 #> 
-#> → calculate ∂(f)/∂(η)
-#> [====|====|====|====|====|====|====|====|====|====] 0:00:01 
+#> Results:
 #> 
-#> → calculate ∂(R²)/∂(η)
-#> [====|====|====|====|====|====|====|====|====|====] 0:00:08 
+#> Likelihood computed by importance sampling
+#>    AIC  BIC logLik
+#>   2276 2272  -1120
 #> 
-#> → finding duplicate expressions in inner model...
-#> [====|====|====|====|====|====|====|====|====|====] 0:00:07 
+#> Optimised parameters:
+#>                 est.   lower   upper
+#> DMTA_0       88.5943 84.3961 92.7925
+#> log_k_DMTA   -3.0466 -3.5609 -2.5322
+#> log_k_M23    -4.0684 -4.9340 -3.2029
+#> log_k_M27    -3.8628 -4.2627 -3.4628
+#> log_k_M31    -3.9803 -4.4804 -3.4801
+#> f_DMTA_ilr_1  0.1304 -0.2186  0.4795
+#> f_DMTA_ilr_2  0.1490 -0.2559  0.5540
+#> f_DMTA_ilr_3 -1.3970 -1.6976 -1.0964
 #> 
-#> → optimizing duplicate expressions in inner model...
-#> [====|====|====|====|====|====|====|====|====|====] 0:00:06 
+#> Correlation: 
+#>              DMTA_0  l__DMTA lg__M23 lg__M27 lg__M31 f_DMTA__1 f_DMTA__2
+#> log_k_DMTA    0.0309                                                    
+#> log_k_M23    -0.0231 -0.0031                                            
+#> log_k_M27    -0.0381 -0.0048  0.0039                                    
+#> log_k_M31    -0.0251 -0.0031  0.0021  0.0830                            
+#> f_DMTA_ilr_1 -0.0046 -0.0006  0.0417 -0.0437  0.0328                    
+#> f_DMTA_ilr_2 -0.0008 -0.0002  0.0214 -0.0270 -0.0909 -0.0361            
+#> f_DMTA_ilr_3 -0.1832 -0.0135  0.0434  0.0804  0.0395 -0.0070    0.0059  
 #> 
-#> → finding duplicate expressions in EBE model...
-#> [====|====|====|====|====|====|====|====|====|====] 0:00:00 
+#> Random effects:
+#>                   est.   lower  upper
+#> SD.DMTA_0       3.3651 -0.9649 7.6951
+#> SD.log_k_DMTA   0.6415  0.2774 1.0055
+#> SD.log_k_M23    1.0176  0.3809 1.6543
+#> SD.log_k_M27    0.4538  0.1522 0.7554
+#> SD.log_k_M31    0.5684  0.1905 0.9464
+#> SD.f_DMTA_ilr_1 0.4111  0.1524 0.6699
+#> SD.f_DMTA_ilr_2 0.4788  0.1808 0.7768
+#> SD.f_DMTA_ilr_3 0.3501  0.1316 0.5685
 #> 
-#> → optimizing duplicate expressions in EBE model...
-#> [====|====|====|====|====|====|====|====|====|====] 0:00:00 
+#> Variance model:
+#>       est.  lower upper
+#> a.1 0.9349 0.8409 1.029
+#> b.1 0.1344 0.1178 0.151
 #> 
-#> → compiling inner model...
-#>  
-#> ✔ done
-#> → finding duplicate expressions in FD model...
+#> Backtransformed parameters:
+#>                   est.     lower    upper
+#> DMTA_0        88.59431 84.396147 92.79246
+#> k_DMTA         0.04752  0.028413  0.07948
+#> k_M23          0.01710  0.007198  0.04064
+#> k_M27          0.02101  0.014084  0.03134
+#> k_M31          0.01868  0.011329  0.03080
+#> f_DMTA_to_M23  0.14498        NA       NA
+#> f_DMTA_to_M27  0.12056        NA       NA
+#> f_DMTA_to_M31  0.11015        NA       NA
 #> 
-#> → optimizing duplicate expressions in FD model...
+#> Resulting formation fractions:
+#>               ff
+#> DMTA_M23  0.1450
+#> DMTA_M27  0.1206
+#> DMTA_M31  0.1101
+#> DMTA_sink 0.6243
 #> 
-#> → compiling EBE model...
-#>  
-#> ✔ done
-#> → compiling events FD model...
-#>  
-#> ✔ done
-#> Model:
-#> cmt(DMTA);
-#> cmt(M23);
-#> cmt(M27);
-#> cmt(M31);
-#> rx_expr_14~ETA[1]+THETA[1];
-#> DMTA(0)=rx_expr_14;
-#> rx_expr_15~ETA[5]+THETA[5];
-#> rx_expr_16~ETA[7]+THETA[7];
-#> rx_expr_17~ETA[6]+THETA[6];
-#> rx_expr_24~exp(rx_expr_15);
-#> rx_expr_25~exp(rx_expr_17);
-#> rx_expr_29~t*rx_expr_24;
-#> rx_expr_30~t*rx_expr_25;
-#> rx_expr_31~exp(-(rx_expr_16));
-#> rx_expr_35~1+rx_expr_31;
-#> rx_expr_40~1/(rx_expr_35);
-#> rx_expr_42~(rx_expr_40);
-#> rx_expr_43~1-rx_expr_42;
-#> d/dt(DMTA)=-DMTA*(exp(rx_expr_15-rx_expr_29)/(rx_expr_35)+exp(rx_expr_17-rx_expr_30)*(rx_expr_43))/(exp(-t*rx_expr_24)/(rx_expr_35)+exp(-t*rx_expr_25)*(rx_expr_43));
-#> rx_expr_18~ETA[2]+THETA[2];
-#> rx_expr_26~exp(rx_expr_18);
-#> d/dt(M23)=-rx_expr_26*M23+DMTA*(exp(rx_expr_15-rx_expr_29)/(rx_expr_35)+exp(rx_expr_17-rx_expr_30)*(rx_expr_43))*f_DMTA_tffm0_1/(exp(-t*rx_expr_24)/(rx_expr_35)+exp(-t*rx_expr_25)*(rx_expr_43));
-#> rx_expr_19~ETA[3]+THETA[3];
-#> rx_expr_20~ETA[4]+THETA[4];
-#> rx_expr_21~1-f_DMTA_tffm0_1;
-#> rx_expr_27~exp(rx_expr_19);
-#> rx_expr_28~exp(rx_expr_20);
-#> d/dt(M27)=-rx_expr_27*M27+rx_expr_28*M31+DMTA*(rx_expr_21)*(exp(rx_expr_15-rx_expr_29)/(rx_expr_35)+exp(rx_expr_17-rx_expr_30)*(rx_expr_43))*f_DMTA_tffm0_2/(exp(-t*rx_expr_24)/(rx_expr_35)+exp(-t*rx_expr_25)*(rx_expr_43));
-#> rx_expr_22~1-f_DMTA_tffm0_2;
-#> d/dt(M31)=-rx_expr_28*M31+DMTA*(rx_expr_22)*(rx_expr_21)*(exp(rx_expr_15-rx_expr_29)/(rx_expr_35)+exp(rx_expr_17-rx_expr_30)*(rx_expr_43))*f_DMTA_tffm0_3/(exp(-t*rx_expr_24)/(rx_expr_35)+exp(-t*rx_expr_25)*(rx_expr_43));
-#> rx_expr_0~CMT==4;
-#> rx_expr_1~CMT==2;
-#> rx_expr_2~CMT==1;
-#> rx_expr_3~CMT==3;
-#> rx_expr_4~1-(rx_expr_0);
-#> rx_expr_5~1-(rx_expr_1);
-#> rx_expr_6~1-(rx_expr_3);
-#> rx_yj_~(rx_expr_4)*((2*(rx_expr_5)*(rx_expr_2)+2*(rx_expr_1))*(rx_expr_6)+2*(rx_expr_3))+2*(rx_expr_0);
-#> rx_expr_7~(rx_expr_1);
-#> rx_expr_8~(rx_expr_3);
-#> rx_expr_9~(rx_expr_0);
-#> rx_expr_13~(rx_expr_5);
-#> rx_expr_32~rx_expr_13*(rx_expr_2);
-#> rx_lambda_~(rx_expr_4)*((rx_expr_32+rx_expr_7)*(rx_expr_6)+rx_expr_8)+rx_expr_9;
-#> rx_hi_~(rx_expr_4)*((rx_expr_32+rx_expr_7)*(rx_expr_6)+rx_expr_8)+rx_expr_9;
-#> rx_low_~0;
-#> rx_expr_10~M31*(rx_expr_0);
-#> rx_expr_11~M27*(rx_expr_3);
-#> rx_expr_12~M23*(rx_expr_1);
-#> rx_expr_23~DMTA*(rx_expr_5);
-#> rx_expr_36~rx_expr_23*(rx_expr_2);
-#> rx_pred_=(rx_expr_4)*((rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))*(rx_expr_3)+((rx_expr_1)*(rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))+(rx_expr_5)*(rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))*(rx_expr_2))*(rx_expr_6))+(rx_expr_0)*(rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)));
-#> rx_expr_33~Rx_pow_di(THETA[12],2);
-#> rx_expr_34~Rx_pow_di(THETA[11],2);
-#> rx_r_=(rx_expr_4)*((rx_expr_33*Rx_pow_di(((rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))*(rx_expr_3)+((rx_expr_1)*(rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))+(rx_expr_5)*(rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))*(rx_expr_2))*(rx_expr_6)),2)+rx_expr_34)*(rx_expr_3)+((rx_expr_1)*(rx_expr_33*Rx_pow_di(((rx_expr_1)*(rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))+(rx_expr_5)*(rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))*(rx_expr_2)),2)+rx_expr_34)+(rx_expr_33*Rx_pow_di(((rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))*(rx_expr_2)),2)+rx_expr_34)*(rx_expr_5)*(rx_expr_2))*(rx_expr_6))+(rx_expr_0)*(rx_expr_33*Rx_pow_di(((rx_expr_4)*((rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))*(rx_expr_3)+((rx_expr_1)*(rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))+(rx_expr_5)*(rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))*(rx_expr_2))*(rx_expr_6))+(rx_expr_0)*(rx_expr_10+(rx_expr_4)*(rx_expr_11+(rx_expr_12+rx_expr_36)*(rx_expr_6)))),2)+rx_expr_34);
-#> DMTA_0=THETA[1];
-#> log_k_M23=THETA[2];
-#> log_k_M27=THETA[3];
-#> log_k_M31=THETA[4];
-#> log_k1=THETA[5];
-#> log_k2=THETA[6];
-#> g_qlogis=THETA[7];
-#> f_DMTA_tffm0_1_qlogis=THETA[8];
-#> f_DMTA_tffm0_2_qlogis=THETA[9];
-#> f_DMTA_tffm0_3_qlogis=THETA[10];
-#> sigma_low=THETA[11];
-#> rsd_high=THETA[12];
-#> eta.DMTA_0=ETA[1];
-#> eta.log_k_M23=ETA[2];
-#> eta.log_k_M27=ETA[3];
-#> eta.log_k_M31=ETA[4];
-#> eta.log_k1=ETA[5];
-#> eta.log_k2=ETA[6];
-#> eta.g_qlogis=ETA[7];
-#> eta.f_DMTA_tffm0_1_qlogis=ETA[8];
-#> eta.f_DMTA_tffm0_2_qlogis=ETA[9];
-#> eta.f_DMTA_tffm0_3_qlogis=ETA[10];
-#> DMTA_0_model=rx_expr_14;
-#> k_M23=rx_expr_26;
-#> k_M27=rx_expr_27;
-#> k_M31=rx_expr_28;
-#> k1=rx_expr_24;
-#> k2=rx_expr_25;
-#> g=1/(rx_expr_35);
-#> f_DMTA_to_M23=1/(1+exp(-(ETA[8]+THETA[8])));
-#> f_DMTA_to_M23=1/(1+exp(-(ETA[9]+THETA[9])));
-#> f_DMTA_to_M23=1/(1+exp(-(ETA[10]+THETA[10])));
-#> f_DMTA_to_M23=f_DMTA_tffm0_1;
-#> f_DMTA_to_M27=(rx_expr_21)*f_DMTA_tffm0_2;
-#> f_DMTA_to_M31=(rx_expr_22)*(rx_expr_21)*f_DMTA_tffm0_3;
-#> tad=tad();
-#> dosenum=dosenum();
-#> Needed Covariates:
-#> [1] "f_DMTA_tffm0_1" "f_DMTA_tffm0_2" "f_DMTA_tffm0_3" "CMT"           
-#> Error in (function (data, inits, PKpars, model = NULL, pred = NULL, err = NULL,     lower = -Inf, upper = Inf, fixed = NULL, skipCov = NULL,     control = foceiControl(), thetaNames = NULL, etaNames = NULL,     etaMat = NULL, ..., env = NULL, keep = NULL, drop = NULL) {    set.seed(control$seed)    .pt <- proc.time()    RxODE::.setWarnIdSort(FALSE)    on.exit(RxODE::.setWarnIdSort(TRUE))    loadNamespace("n1qn1")    if (!RxODE::rxIs(control, "foceiControl")) {        control <- do.call(foceiControl, control)    }    if (is.null(env)) {        .ret <- new.env(parent = emptyenv())    }    else {        .ret <- env    }    .ret$origData <- data    .ret$etaNames <- etaNames    .ret$thetaFixed <- fixed    .ret$control <- control    .ret$control$focei.mu.ref <- integer(0)    if (is(model, "RxODE") || is(model, "character")) {        .ret$ODEmodel <- TRUE        if (class(pred) != "function") {            stop("pred must be a function specifying the prediction variables in this model.")        }    }    else {        .ret$ODEmodel <- TRUE        model <- RxODE::rxGetLin(PKpars)        pred <- eval(parse(text = "function(){return(Central);}"))    }    .square <- function(x) x * x    .ret$diagXformInv <- c(sqrt = ".square", log = "exp", identity = "identity")[control$diagXform]    if (is.null(err)) {        err <- eval(parse(text = paste0("function(){err", paste(inits$ERROR[[1]],             collapse = ""), "}")))    }    .covNames <- .parNames <- c()    .ret$adjLik <- control$adjLik    .mixed <- !is.null(inits$OMGA) && length(inits$OMGA) > 0    if (!exists("noLik", envir = .ret)) {        .atol <- rep(control$atol, length(RxODE::rxModelVars(model)$state))        .rtol <- rep(control$rtol, length(RxODE::rxModelVars(model)$state))        .ssAtol <- rep(control$ssAtol, length(RxODE::rxModelVars(model)$state))        .ssRtol <- rep(control$ssRtol, length(RxODE::rxModelVars(model)$state))        .ret$model <- RxODE::rxSymPySetupPred(model, pred, PKpars,             err, grad = (control$derivMethod == 2L), pred.minus.dv = TRUE,             sum.prod = control$sumProd, theta.derivs = FALSE,             optExpression = control$optExpression, interaction = (control$interaction ==                 1L), only.numeric = !.mixed, run.internal = TRUE,             addProp = control$addProp)        if (!is.null(.ret$model$inner)) {            .atol <- c(.atol, rep(control$atolSens, length(RxODE::rxModelVars(.ret$model$inner)$state) -                 length(.atol)))            .rtol <- c(.rtol, rep(control$rtolSens, length(RxODE::rxModelVars(.ret$model$inner)$state) -                 length(.rtol)))            .ret$control$rxControl$atol <- .atol            .ret$control$rxControl$rtol <- .rtol            .ssAtol <- c(.ssAtol, rep(control$ssAtolSens, length(RxODE::rxModelVars(.ret$model$inner)$state) -                 length(.ssAtol)))            .ssRtol <- c(.ssRtol, rep(control$ssRtolSens, length(RxODE::rxModelVars(.ret$model$inner)$state) -                 length(.ssRtol)))            .ret$control$rxControl$ssAtol <- .ssAtol            .ret$control$rxControl$ssRtol <- .ssRtol        }        .covNames <- .parNames <- RxODE::rxParams(.ret$model$pred.only)        .covNames <- .covNames[regexpr(rex::rex(start, or("THETA",             "ETA"), "[", numbers, "]", end), .covNames) == -1]        colnames(data) <- sapply(names(data), function(x) {            if (any(x == .covNames)) {                return(x)            }            else {                return(toupper(x))            }        })        .lhs <- c(names(RxODE::rxInits(.ret$model$pred.only)),             RxODE::rxLhs(.ret$model$pred.only))        if (length(.lhs) > 0) {            .covNames <- .covNames[regexpr(rex::rex(start, or(.lhs),                 end), .covNames) == -1]        }        if (length(.covNames) > 0) {            if (!all(.covNames %in% names(data))) {                message("Model:")                RxODE::rxCat(.ret$model$pred.only)                message("Needed Covariates:")                nlmixrPrint(.covNames)                stop("Not all the covariates are in the dataset.")            }            message("Needed Covariates:")            print(.covNames)        }        .extraPars <- .ret$model$extra.pars    }    else {        if (.ret$noLik) {            .atol <- rep(control$atol, length(RxODE::rxModelVars(model)$state))            .rtol <- rep(control$rtol, length(RxODE::rxModelVars(model)$state))            .ret$model <- RxODE::rxSymPySetupPred(model, pred,                 PKpars, err, grad = FALSE, pred.minus.dv = TRUE,                 sum.prod = control$sumProd, theta.derivs = FALSE,                 optExpression = control$optExpression, run.internal = TRUE,                 only.numeric = TRUE, addProp = control$addProp)            if (!is.null(.ret$model$inner)) {                .atol <- c(.atol, rep(control$atolSens, length(RxODE::rxModelVars(.ret$model$inner)$state) -                   length(.atol)))                .rtol <- c(.rtol, rep(control$rtolSens, length(RxODE::rxModelVars(.ret$model$inner)$state) -                   length(.rtol)))                .ret$control$rxControl$atol <- .atol                .ret$control$rxControl$rtol <- .rtol            }            .covNames <- .parNames <- RxODE::rxParams(.ret$model$pred.only)            .covNames <- .covNames[regexpr(rex::rex(start, or("THETA",                 "ETA"), "[", numbers, "]", end), .covNames) ==                 -1]            colnames(data) <- sapply(names(data), function(x) {                if (any(x == .covNames)) {                  return(x)                }                else {                  return(toupper(x))                }            })            .lhs <- c(names(RxODE::rxInits(.ret$model$pred.only)),                 RxODE::rxLhs(.ret$model$pred.only))            if (length(.lhs) > 0) {                .covNames <- .covNames[regexpr(rex::rex(start,                   or(.lhs), end), .covNames) == -1]            }            if (length(.covNames) > 0) {                if (!all(.covNames %in% names(data))) {                  message("Model:")                  RxODE::rxCat(.ret$model$pred.only)                  message("Needed Covariates:")                  nlmixrPrint(.covNames)                  stop("Not all the covariates are in the dataset.")                }                message("Needed Covariates:")                print(.covNames)            }            .extraPars <- .ret$model$extra.pars        }        else {            .extraPars <- NULL        }    }    .ret$skipCov <- skipCov    if (is.null(skipCov)) {        if (is.null(fixed)) {            .tmp <- rep(FALSE, length(inits$THTA))        }        else {            if (length(fixed) < length(inits$THTA)) {                .tmp <- c(fixed, rep(FALSE, length(inits$THTA) -                   length(fixed)))            }            else {                .tmp <- fixed[1:length(inits$THTA)]            }        }        if (exists("uif", envir = .ret)) {            .uifErr <- .ret$uif$ini$err[!is.na(.ret$uif$ini$ntheta)]            .uifErr <- sapply(.uifErr, function(x) {                if (is.na(x)) {                  return(FALSE)                }                return(!any(x == c("pow2", "tbs", "tbsYj")))            })            .tmp <- (.tmp | .uifErr)        }        .ret$skipCov <- c(.tmp, rep(TRUE, length(.extraPars)))        .ret$control$focei.mu.ref <- .ret$uif$focei.mu.ref    }    if (is.null(.extraPars)) {        .nms <- c(sprintf("THETA[%s]", seq_along(inits$THTA)))    }    else {        .nms <- c(sprintf("THETA[%s]", seq_along(inits$THTA)),             sprintf("ERR[%s]", seq_along(.extraPars)))    }    if (!is.null(thetaNames) && (length(inits$THTA) + length(.extraPars)) ==         length(thetaNames)) {        .nms <- thetaNames    }    .ret$thetaNames <- .nms    .thetaReset$thetaNames <- .nms    if (length(lower) == 1) {        lower <- rep(lower, length(inits$THTA))    }    else if (length(lower) != length(inits$THTA)) {        print(inits$THTA)        print(lower)        stop("Lower must be a single constant for all the THETA lower bounds, or match the dimension of THETA.")    }    if (length(upper) == 1) {        upper <- rep(upper, length(inits$THTA))    }    else if (length(lower) != length(inits$THTA)) {        stop("Upper must be a single constant for all the THETA lower bounds, or match the dimension of THETA.")    }    if (!is.null(.extraPars)) {        .ret$model$extra.pars <- eval(call(control$diagXform,             .ret$model$extra.pars))        if (length(.ret$model$extra.pars) > 0) {            inits$THTA <- c(inits$THTA, .ret$model$extra.pars)            .lowerErr <- rep(control$atol[1] * 10, length(.ret$model$extra.pars))            .upperErr <- rep(Inf, length(.ret$model$extra.pars))            lower <- c(lower, .lowerErr)            upper <- c(upper, .upperErr)        }    }    if (is.null(data$ID))         stop("\"ID\" not found in data")    if (is.null(data$DV))         stop("\"DV\" not found in data")    if (is.null(data$EVID))         data$EVID <- 0    if (is.null(data$AMT))         data$AMT <- 0    for (.v in c("TIME", "AMT", "DV", .covNames)) {        data[[.v]] <- as.double(data[[.v]])    }    .ret$dataSav <- data    .ds <- data[data$EVID != 0 & data$EVID != 2, c("ID", "TIME",         "AMT", "EVID", .covNames)]    .w <- which(tolower(names(data)) == "limit")    .limitName <- NULL    if (length(.w) == 1L) {        .limitName <- names(data)[.w]    }    .censName <- NULL    .w <- which(tolower(names(data)) == "cens")    if (length(.w) == 1L) {        .censName <- names(data[.w])    }    data <- data[data$EVID == 0 | data$EVID == 2, c("ID", "TIME",         "DV", "EVID", .covNames, .limitName, .censName)]    .w <- which(!(names(.ret$dataSav) %in% c(.covNames, keep)))    names(.ret$dataSav)[.w] <- tolower(names(.ret$dataSav[.w]))    if (.mixed) {        .lh <- .parseOM(inits$OMGA)        .nlh <- sapply(.lh, length)        .osplt <- rep(1:length(.lh), .nlh)        .lini <- list(inits$THTA, unlist(.lh))        .nlini <- sapply(.lini, length)        .nsplt <- rep(1:length(.lini), .nlini)        .om0 <- .genOM(.lh)        if (length(etaNames) == dim(.om0)[1]) {            .ret$etaNames <- .ret$etaNames        }        else {            .ret$etaNames <- sprintf("ETA[%d]", seq(1, dim(.om0)[1]))        }        .ret$rxInv <- RxODE::rxSymInvCholCreate(mat = .om0, diag.xform = control$diagXform)        .ret$xType <- .ret$rxInv$xType        .om0a <- .om0        .om0a <- .om0a/control$diagOmegaBoundLower        .om0b <- .om0        .om0b <- .om0b * control$diagOmegaBoundUpper        .om0a <- RxODE::rxSymInvCholCreate(mat = .om0a, diag.xform = control$diagXform)        .om0b <- RxODE::rxSymInvCholCreate(mat = .om0b, diag.xform = control$diagXform)        .omdf <- data.frame(a = .om0a$theta, m = .ret$rxInv$theta,             b = .om0b$theta, diag = .om0a$theta.diag)        .omdf$lower <- with(.omdf, ifelse(a > b, b, a))        .omdf$lower <- with(.omdf, ifelse(lower == m, -Inf, lower))        .omdf$lower <- with(.omdf, ifelse(!diag, -Inf, lower))        .omdf$upper <- with(.omdf, ifelse(a < b, b, a))        .omdf$upper <- with(.omdf, ifelse(upper == m, Inf, upper))        .omdf$upper <- with(.omdf, ifelse(!diag, Inf, upper))        .ret$control$nomega <- length(.omdf$lower)        .ret$control$neta <- sum(.omdf$diag)        .ret$control$ntheta <- length(lower)        .ret$control$nfixed <- sum(fixed)        lower <- c(lower, .omdf$lower)        upper <- c(upper, .omdf$upper)    }    else {        .ret$control$nomega <- 0        .ret$control$neta <- 0        .ret$xType <- -1        .ret$control$ntheta <- length(lower)        .ret$control$nfixed <- sum(fixed)    }    .ret$lower <- lower    .ret$upper <- upper    .ret$thetaIni <- inits$THTA    .scaleC <- double(length(lower))    if (is.null(control$scaleC)) {        .scaleC <- rep(NA_real_, length(lower))    }    else {        .scaleC <- as.double(control$scaleC)        if (length(lower) > length(.scaleC)) {            .scaleC <- c(.scaleC, rep(NA_real_, length(lower) -                 length(.scaleC)))        }        else if (length(lower) < length(.scaleC)) {            .scaleC <- .scaleC[seq(1, length(lower))]            warning("scaleC control option has more options than estimated population parameters, please check.")        }    }    .ret$scaleC <- .scaleC    if (exists("uif", envir = .ret)) {        .ini <- as.data.frame(.ret$uif$ini)[!is.na(.ret$uif$ini$err),             c("est", "err", "ntheta")]        for (.i in seq_along(.ini$err)) {            if (is.na(.ret$scaleC[.ini$ntheta[.i]])) {                if (any(.ini$err[.i] == c("boxCox", "yeoJohnson",                   "pow2", "tbs", "tbsYj"))) {                  .ret$scaleC[.ini$ntheta[.i]] <- 1                }                else if (any(.ini$err[.i] == c("prop", "add",                   "norm", "dnorm", "logn", "dlogn", "lnorm",                   "dlnorm"))) {                  .ret$scaleC[.ini$ntheta[.i]] <- 0.5 * abs(.ini$est[.i])                }            }        }        for (.i in .ini$model$extraProps$powTheta) {            if (is.na(.ret$scaleC[.i]))                 .ret$scaleC[.i] <- 1        }        .ini <- as.data.frame(.ret$uif$ini)        for (.i in .ini$model$extraProps$factorial) {            if (is.na(.ret$scaleC[.i]))                 .ret$scaleC[.i] <- abs(1/digamma(.ini$est[.i] +                   1))        }        for (.i in .ini$model$extraProps$gamma) {            if (is.na(.ret$scaleC[.i]))                 .ret$scaleC[.i] <- abs(1/digamma(.ini$est[.i]))        }        for (.i in .ini$model$extraProps$log) {            if (is.na(.ret$scaleC[.i]))                 .ret$scaleC[.i] <- log(abs(.ini$est[.i])) * abs(.ini$est[.i])        }        for (.i in .ret$logitThetas) {            .b <- .ret$logitThetasLow[.i]            .c <- .ret$logitThetasHi[.i]            .a <- .ini$est[.i]            if (is.na(.ret$scaleC[.i])) {                .ret$scaleC[.i] <- 1 * (-.b + .c) * exp(-.a)/((1 +                   exp(-.a))^2 * (.b + 1 * (-.b + .c)/(1 + exp(-.a))))            }        }    }    names(.ret$thetaIni) <- sprintf("THETA[%d]", seq_along(.ret$thetaIni))    if (is.null(etaMat) & !is.null(control$etaMat)) {        .ret$etaMat <- control$etaMat    }    else {        .ret$etaMat <- etaMat    }    .ret$setupTime <- (proc.time() - .pt)["elapsed"]    if (exists("uif", envir = .ret)) {        .tmp <- .ret$uif$logThetasList        .ret$logThetas <- .tmp[[1]]        .ret$logThetasF <- .tmp[[2]]        .tmp <- .ret$uif$logitThetasList        .ret$logitThetas <- .tmp[[1]]        .ret$logitThetasF <- .tmp[[2]]        .tmp <- .ret$uif$logitThetasListLow        .ret$logitThetasLow <- .tmp[[1]]        .ret$logitThetasLowF <- .tmp[[2]]        .tmp <- .ret$uif$logitThetasListHi        .ret$logitThetasHi <- .tmp[[1]]        .ret$logitThetasHiF <- .tmp[[2]]        .tmp <- .ret$uif$probitThetasList        .ret$probitThetas <- .tmp[[1]]        .ret$probitThetasF <- .tmp[[2]]        .tmp <- .ret$uif$probitThetasListLow        .ret$probitThetasLow <- .tmp[[1]]        .ret$probitThetasLowF <- .tmp[[2]]        .tmp <- .ret$uif$probitThetasListHi        .ret$probitThetasHi <- .tmp[[1]]        .ret$probitThetasHiF <- .tmp[[2]]    }    else {        .ret$logThetasF <- integer(0)        .ret$logitThetasF <- integer(0)        .ret$logitThetasHiF <- numeric(0)        .ret$logitThetasLowF <- numeric(0)        .ret$logitThetas <- integer(0)        .ret$logitThetasHi <- numeric(0)        .ret$logitThetasLow <- numeric(0)        .ret$probitThetasF <- integer(0)        .ret$probitThetasHiF <- numeric(0)        .ret$probitThetasLowF <- numeric(0)        .ret$probitThetas <- integer(0)        .ret$probitThetasHi <- numeric(0)        .ret$probitThetasLow <- numeric(0)    }    if (exists("noLik", envir = .ret)) {        if (!.ret$noLik) {            .ret$.params <- c(sprintf("THETA[%d]", seq_along(.ret$thetaIni)),                 sprintf("ETA[%d]", seq(1, dim(.om0)[1])))            .ret$.thetan <- length(.ret$thetaIni)            .ret$nobs <- sum(data$EVID == 0)        }    }    .ret$control$printTop <- TRUE    .ret$control$nF <- 0    .est0 <- .ret$thetaIni    if (!is.null(.ret$model$pred.nolhs)) {        .ret$control$predNeq <- length(.ret$model$pred.nolhs$state)    }    else {        .ret$control$predNeq <- 0L    }    .fitFun <- function(.ret) {        this.env <- environment()        assign("err", "theta reset", this.env)        while (this.env$err == "theta reset") {            assign("err", "", this.env)            .ret0 <- tryCatch({                foceiFitCpp_(.ret)            }, error = function(e) {                if (regexpr("theta reset", e$message) != -1) {                  assign("zeroOuter", FALSE, this.env)                  assign("zeroGrad", FALSE, this.env)                  if (regexpr("theta reset0", e$message) != -1) {                    assign("zeroGrad", TRUE, this.env)                  }                  else if (regexpr("theta resetZ", e$message) !=                     -1) {                    assign("zeroOuter", TRUE, this.env)                  }                  assign("err", "theta reset", this.env)                }                else {                  assign("err", e$message, this.env)                }            })            if (this.env$err == "theta reset") {                .nm <- names(.ret$thetaIni)                .ret$thetaIni <- setNames(.thetaReset$thetaIni +                   0, .nm)                .ret$rxInv$theta <- .thetaReset$omegaTheta                .ret$control$printTop <- FALSE                .ret$etaMat <- .thetaReset$etaMat                .ret$control$etaMat <- .thetaReset$etaMat                .ret$control$maxInnerIterations <- .thetaReset$maxInnerIterations                .ret$control$nF <- .thetaReset$nF                .ret$control$gillRetC <- .thetaReset$gillRetC                .ret$control$gillRet <- .thetaReset$gillRet                .ret$control$gillRet <- .thetaReset$gillRet                .ret$control$gillDf <- .thetaReset$gillDf                .ret$control$gillDf2 <- .thetaReset$gillDf2                .ret$control$gillErr <- .thetaReset$gillErr                .ret$control$rEps <- .thetaReset$rEps                .ret$control$aEps <- .thetaReset$aEps                .ret$control$rEpsC <- .thetaReset$rEpsC                .ret$control$aEpsC <- .thetaReset$aEpsC                .ret$control$c1 <- .thetaReset$c1                .ret$control$c2 <- .thetaReset$c2                if (this.env$zeroOuter) {                  message("Posthoc reset")                  .ret$control$maxOuterIterations <- 0L                }                else if (this.env$zeroGrad) {                  message("Theta reset (zero gradient values); Switch to bobyqa")                  RxODE::rxReq("minqa")                  .ret$control$outerOptFun <- .bobyqa                  .ret$control$outerOpt <- -1L                }                else {                  message("Theta reset (ETA drift)")                }            }        }        if (this.env$err != "") {            stop(this.env$err)        }        else {            return(.ret0)        }    }    .ret0 <- try(.fitFun(.ret))    .n <- 1    while (inherits(.ret0, "try-error") && control$maxOuterIterations !=         0 && .n <= control$nRetries) {        message(sprintf("Restart %s", .n))        .ret$control$nF <- 0        .estNew <- .est0 + 0.2 * .n * abs(.est0) * stats::runif(length(.est0)) -             0.1 * .n        .estNew <- sapply(seq_along(.est0), function(.i) {            if (.ret$thetaFixed[.i]) {                return(.est0[.i])            }            else if (.estNew[.i] < lower[.i]) {                return(lower + (.Machine$double.eps)^(1/7))            }            else if (.estNew[.i] > upper[.i]) {                return(upper - (.Machine$double.eps)^(1/7))            }            else {                return(.estNew[.i])            }        })        .ret$thetaIni <- .estNew        .ret0 <- try(.fitFun(.ret))        .n <- .n + 1    }    if (inherits(.ret0, "try-error"))         stop("Could not fit data.")    .ret <- .ret0    if (exists("parHistData", .ret)) {        .tmp <- .ret$parHistData        .tmp <- .tmp[.tmp$type == "Unscaled", names(.tmp) !=             "type"]        .iter <- .tmp$iter        .tmp <- .tmp[, names(.tmp) != "iter"]        .ret$parHistStacked <- data.frame(stack(.tmp), iter = .iter)        names(.ret$parHistStacked) <- c("val", "par", "iter")        .ret$parHist <- data.frame(iter = .iter, .tmp)    }    if (.mixed) {        .etas <- .ret$ranef        .thetas <- .ret$fixef        .pars <- .Call(`_nlmixr_nlmixrParameters`, .thetas, .etas)        .ret$shrink <- .Call(`_nlmixr_calcShrinkOnly`, .ret$omega,             .pars$eta.lst, length(.etas$ID))        .updateParFixed(.ret)    }    else {        .updateParFixed(.ret)    }    if (!exists("table", .ret)) {        .ret$table <- tableControl()    }    if (control$calcTables) {        .ret <- addTable(.ret, updateObject = "no", keep = keep,             drop = drop, table = .ret$table)    }    .ret})(data = dat, inits = .FoceiInits, PKpars = .pars, model = .mod,     pred = function() {        return(nlmixr_pred)    }, err = uif$error, lower = uif$focei.lower, upper = uif$focei.upper,     fixed = uif$focei.fixed, thetaNames = uif$focei.names, etaNames = uif$eta.names,     control = control, env = env, keep = .keep, drop = .drop): Not all the covariates are in the dataset.
-#> Timing stopped at: 119.8 9.331 129.2
-#> Timing stopped at: 120 9.331 129.3
-summary(f_dmta_nlmixr_focei)
-#> Error in summary(f_dmta_nlmixr_focei): object 'f_dmta_nlmixr_focei' not found
-plot(f_dmta_nlmixr_focei)
-#> Error in plot(f_dmta_nlmixr_focei): object 'f_dmta_nlmixr_focei' not found
-# Using saemix takes about 18 minutes
-system.time(
-  f_dmta_saemix <- saem(f_dmta_mkin_tc, test_log_parms = TRUE)
-)
-#> DINTDY-  T (=R1) illegal      
-#> In above message, R1 = 115.507
-#>  
-#>       T not in interval TCUR - HU (= R1) to TCUR (=R2)      
-#> In above message, R1 = 112.133, R2 = 113.577
-#>  
-#> DLSODA-  At T (=R1), too much accuracy requested  
-#>       for precision of machine..  See TOLSF (=R2) 
-#> In above message, R1 = 55.3899, R2 = nan
-#>  
-#> Error in out[available, var]: (subscript) logical subscript too long
-#> Timing stopped at: 11.84 0.008 11.85
-#> Timing stopped at: 12.16 0.008 12.17
-
-# nlmixr with est = "saem" is pretty fast with default iteration numbers, most
-# of the time (about 2.5 minutes) is spent for calculating the log likelihood at the end
-# The likelihood calculated for the nlmixr fit is much lower than that found by saemix
-# Also, the trace plot and the plot of the individual predictions is not
-# convincing for the parent. It seems we are fitting an overparameterised
-# model, so the result we get strongly depends on starting parameters and control settings.
-system.time(
-  f_dmta_nlmixr_saem <- nlmixr(f_dmta_mkin_tc, est = "saem",
-    control = nlmixr::saemControl(print = 500, logLik = TRUE, nmc = 9))
-)
-#> With est = 'saem', a different error model is required for each observed variableChanging the error model to 'obs_tc' (Two-component error for each observed variable)
-#> ℹ parameter labels from comments are typically ignored in non-interactive mode
-#> ℹ Need to run with the source intact to parse comments
-#> Error in eval(substitute(expr), data, enclos = parent.frame()): Cannot run SAEM since some of the parameters are not mu-referenced (eta.f_DMTA_tffm0_1_qlogis, eta.f_DMTA_tffm0_2_qlogis, eta.f_DMTA_tffm0_3_qlogis)
-#> Timing stopped at: 0.892 0.004 0.896
-#> Timing stopped at: 1.096 0.005 1.1
-traceplot(f_dmta_nlmixr_saem$nm)
-#> Error in traceplot(f_dmta_nlmixr_saem$nm): could not find function "traceplot"
-summary(f_dmta_nlmixr_saem)
-#> Error in summary(f_dmta_nlmixr_saem): object 'f_dmta_nlmixr_saem' not found
-plot(f_dmta_nlmixr_saem)
-#> Error in plot(f_dmta_nlmixr_saem): object 'f_dmta_nlmixr_saem' not found
-# }
+#> Estimated disappearance times:
+#>       DT50   DT90
+#> DMTA 14.59  48.45
+#> M23  40.52 134.62
+#> M27  32.99 109.60
+#> M31  37.11 123.26
+# As the confidence interval for the random effects of DMTA_0
+# includes zero, we could try an alternative model without
+# such random effects
+# f_dmta_saem_tc_2 <- saem(dmta_sfo_sfo3p_tc,
+#   covariance.model = diag(c(0, rep(1, 7))))
+# saemix::plot(f_dmta_saem_tc_2$so, plot.type = "convergence")
+# This does not perform better judged by AIC and BIC
+# saemix::compare.saemix(f_dmta_saem_tc$so, f_dmta_saem_tc_2$so)
+# }
    @@ -446,7 +321,7 @@ specific pieces of information in the comments.

    -

    Site built with pkgdown 2.0.2.

    +

    Site built with pkgdown 2.0.6.

    diff --git a/docs/dev/reference/endpoints.html b/docs/dev/reference/endpoints.html index a472b054..ed72ec47 100644 --- a/docs/dev/reference/endpoints.html +++ b/docs/dev/reference/endpoints.html @@ -23,7 +23,7 @@ advantage that the SFORB model can also be used for metabolites."> mkin - 1.1.0 + 1.1.2 @@ -32,7 +32,7 @@ advantage that the SFORB model can also be used for metabolites.">Functions and data
  • -
    • +
  • + Example evaluations of dimethenamid data from 2018 with nonlinear mixed-effects models +
  • Example evaluation of FOCUS Example Dataset Z
    • @@ -95,21 +98,24 @@ advantage that the SFORB model can also be used for metabolites.

    - 
    endpoints(fit)
    + 
    endpoints(fit)

    Arguments

    fit
    -

    An object of class mkinfit, nlme.mmkin, saem.mmkin or -nlmixr.mmkin. Or another object that has list components +

    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.

    +

    Value

    -

    A list with a matrix of dissipation times named distimes, + + +

    A list with a matrix of dissipation times named distimes, 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

    @@ -131,22 +137,22 @@ HS and DFOP, as well as from Eigenvalues b1 and b2 of any SFORB models

    Examples

    - 
    
-  fit <- mkinfit("FOMC", FOCUS_2006_C, quiet = TRUE)
-  endpoints(fit)
+    
    
+  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)
+  # \dontrun{
+    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)
+    fit_3 <- mkinfit("SFORB", FOCUS_2006_C, quiet = TRUE)
+    endpoints(fit_3)
 #> $ff
 #> parent_free 
 #>           1 
@@ -159,8 +165,8 @@ HS and DFOP, as well as from Eigenvalues b1 and b2 of any SFORB models
    
 #>            DT50     DT90 DT50back DT50_parent_b1 DT50_parent_b2
 #> parent 1.886925 21.25106 6.397208       1.508293       38.83438
 #> 
-  # }
-
+  # }
+
    @@ -175,7 +181,7 @@ HS and DFOP, as well as from Eigenvalues b1 and b2 of any SFORB models

    -

    Site built with pkgdown 2.0.2.

    +

    Site built with pkgdown 2.0.6.

    diff --git a/docs/dev/reference/illparms.html b/docs/dev/reference/illparms.html new file mode 100644 index 00000000..184299ed --- /dev/null +++ b/docs/dev/reference/illparms.html @@ -0,0 +1,207 @@ + +Method to get the names of ill-defined parameters — illparms • mkin + + +
    +
    + + + +
    +
    +
    +

    Method to get the names of ill-defined parameters

    + Source: R/illparms.R + +
    + +
    +

    The method for generalised nonlinear regression fits as obtained +with mkinfit and mmkin checks if the degradation parameters +pass the Wald test (in degradation kinetics often simply called t-test) for +significant difference from zero. For this test, the parameterisation +without parameter transformations is used.

    +
    + +
    + 
    illparms(object, ...)
+
+# S3 method for mkinfit
+illparms(object, conf.level = 0.95, ...)
+
+# S3 method for mmkin
+illparms(object, conf.level = 0.95, ...)
+
+# S3 method for illparms.mmkin
+print(x, ...)
+
+# S3 method for saem.mmkin
+illparms(object, conf.level = 0.95, random = TRUE, errmod = TRUE, ...)
+
+# S3 method for mhmkin
+illparms(object, conf.level = 0.95, random = TRUE, errmod = TRUE, ...)
+
+# S3 method for illparms.mhmkin
+print(x, ...)
    +
    + +
    +

    Arguments

    +
    object
    +

    The object to investigate

    + + +
    ...
    +

    For potential future extensions

    + + +
    conf.level
    +

    The confidence level for checking p values

    + + +
    x
    +

    The object to be printed

    + + +
    random
    +

    For hierarchical fits, should random effects be tested?

    + + +
    errmod
    +

    For hierarchical fits, should error model parameters be +tested?

    + +
    +
    +

    Value

    + + +

    For mkinfit or saem objects, a character vector of parameter +names. For mmkin or mhmkin objects, a matrix like object of class +'illparms.mmkin' or 'illparms.mhmkin'. The latter objects have a suitable +printing method.

    +
    +
    +

    Details

    +

    The method for hierarchical model fits, also known as nonlinear +mixed-effects model fits as obtained with saem and mhmkin +checks if any of the confidence intervals for the random +effects expressed as standard deviations include zero, and if +the confidence intervals for the error model parameters include +zero.

    +
    + +
    +

    Examples

    + 
    fit <- mkinfit("FOMC", FOCUS_2006_A, quiet = TRUE)
+illparms(fit)
+#> [1] "parent_0" "alpha"    "beta"     "sigma"   
+# \dontrun{
+fits <- mmkin(
+  c("SFO", "FOMC"),
+  list("FOCUS A" = FOCUS_2006_A,
+       "FOCUS C" = FOCUS_2006_C),
+  quiet = TRUE)
+illparms(fits)
+#>       dataset
+#> model  FOCUS A                      FOCUS C
+#>   SFO                                      
+#>   FOMC parent_0, alpha, beta, sigma        
+# }
+
    +
    +
    + +
    + + +
    +

    Developed by Johannes Ranke.

    +
    + +
    +

    Site built with pkgdown 2.0.6.

    +
    + +
    + + + + + + + + diff --git a/docs/dev/reference/index.html b/docs/dev/reference/index.html index afd3692d..98f170b6 100644 --- a/docs/dev/reference/index.html +++ b/docs/dev/reference/index.html @@ -17,7 +17,7 @@ mkin - 1.1.0 + 1.1.2 @@ -26,7 +26,7 @@ Functions and data
  • -
  • Example evaluation of FOCUS Laboratory Data L1 to L3
    • +
  • + Example evaluations of dimethenamid data from 2018 with nonlinear mixed-effects models +
  • Example evaluation of FOCUS Example Dataset Z
    • @@ -95,9 +98,27 @@

    Fit one or more kinetic models with one or more state variables to one or more datasets

    -

    nafta() print(<nafta>)

    +

    mhmkin() `[`(<mhmkin>) print(<mhmkin>)

    -

    Evaluate parent kinetics using the NAFTA guidance

    +

    Fit nonlinear mixed-effects models built from one or more kinetic +degradation models and one or more error models

    + +

    Generics

    +

    Generic functions introduced by the package

    + + +

    convergence() print(<convergence.mmkin>)

    + +

    Method to get convergence information

    + +

    illparms() print(<illparms.mmkin>) print(<illparms.mhmkin>)

    + +

    Method to get the names of ill-defined parameters

    + +

    endpoints()

    + +

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

    Show results

    Functions working with mkinfit objects

    @@ -134,11 +155,6 @@ more datasets

    mkinerrmin()

    Calculate the minimum error to assume in order to pass the variance test

    - -

    endpoints()

    - -

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

    aw()

    @@ -164,9 +180,13 @@ of an mmkin object

    AIC(<mmkin>) BIC(<mmkin>)

    Calculate the AIC for a column of an mmkin object

    + +

    summary(<mmkin>) print(<summary.mmkin>)

    + +

    Summary method for class "mmkin"

    Mixed models

    -

    Create and work with nonlinear mixed effects models

    +

    Create and work with nonlinear hierarchical models

    nlme(<mmkin>) print(<nlme.mmkin>) update(<nlme.mmkin>)

    @@ -177,9 +197,10 @@ of an mmkin object

    Fit nonlinear mixed models with SAEM

    -

    nlmixr(<mmkin>) print(<nlmixr.mmkin>) nlmixr_model() nlmixr_data()

    +

    mhmkin() `[`(<mhmkin>) print(<mhmkin>)

    -

    Fit nonlinear mixed models using nlmixr

    +

    Fit nonlinear mixed-effects models built from one or more kinetic +degradation models and one or more error models

    plot(<mixed.mmkin>)

    @@ -188,10 +209,6 @@ of an mmkin object

    summary(<nlme.mmkin>) print(<summary.nlme.mmkin>)

    Summary method for class "nlme.mmkin"

    - -

    summary(<nlmixr.mmkin>) print(<summary.nlmixr.mmkin>)

    - -

    Summary method for class "nlmixr.mmkin"

    summary(<saem.mmkin>) print(<summary.saem.mmkin>)

    @@ -209,17 +226,13 @@ of an mmkin object

    Create a mixed effects model from an mmkin row object

    -

    reexports

    +

    reexports intervals lrtest nlme

    Objects exported from other packages

    intervals(<saem.mmkin>)

    Confidence intervals for parameters in saem.mmkin objects

    - -

    intervals(<nlmixr.mmkin>)

    - -

    Confidence intervals for parameters in nlmixr.mmkin objects

    Datasets and known results

    @@ -341,10 +354,6 @@ kinetic models fitted with mkinfit

    ilr() invilr()

    Function to perform isometric log-ratio transformation

    - -

    tffm0() invtffm0()

    - -

    Transform formation fractions as in the first published mkin version

    logLik(<mkinfit>)

    @@ -442,7 +451,7 @@ kinetic models fitted with mkinfit

    -

    Site built with pkgdown 2.0.2.

    +

    Site built with pkgdown 2.0.6.

    diff --git a/docs/dev/reference/intervals.saem.mmkin.html b/docs/dev/reference/intervals.saem.mmkin.html index 8e6a348b..ee714ad0 100644 --- a/docs/dev/reference/intervals.saem.mmkin.html +++ b/docs/dev/reference/intervals.saem.mmkin.html @@ -17,7 +17,7 @@ mkin - 1.1.0 + 1.1.2 @@ -26,7 +26,7 @@ Functions and data
  • -
  • Example evaluation of FOCUS Laboratory Data L1 to L3
    • +
  • + Example evaluations of dimethenamid data from 2018 with nonlinear mixed-effects models +
  • Example evaluation of FOCUS Example Dataset Z
    • @@ -84,27 +87,36 @@
    - 
    # S3 method for saem.mmkin
-intervals(object, level = 0.95, backtransform = TRUE, ...)
    + 
    # S3 method for saem.mmkin
+intervals(object, level = 0.95, backtransform = TRUE, ...)

    Arguments

    object

    The fitted saem.mmkin object

    + +
    level

    The confidence level. Must be the default of 0.95 as this is what is available in the saemix object

    + +
    backtransform

    In case the model was fitted with mkin transformations, should we backtransform the parameters where a one to one correlation between transformed and backtransformed parameters exists?

    + +
    ...

    For compatibility with the generic method

    +

    Value

    -

    An object with 'intervals.saem.mmkin' and 'intervals.lme' in the + + +

    An object with 'intervals.saem.mmkin' and 'intervals.lme' in the class attribute

    @@ -120,7 +132,7 @@ class attribute

    -

    Site built with pkgdown 2.0.2.

    +

    Site built with pkgdown 2.0.6.

    diff --git a/docs/dev/reference/mean_degparms.html b/docs/dev/reference/mean_degparms.html index d348fa74..67db1868 100644 --- a/docs/dev/reference/mean_degparms.html +++ b/docs/dev/reference/mean_degparms.html @@ -17,7 +17,7 @@ mkin - 1.1.0 + 1.1.2 @@ -26,7 +26,7 @@ Functions and data
  • -
  • Example evaluation of FOCUS Laboratory Data L1 to L3
    • +
  • + Example evaluations of dimethenamid data from 2018 with nonlinear mixed-effects models +
  • Example evaluation of FOCUS Example Dataset Z
    • @@ -84,26 +87,47 @@
    - 
    mean_degparms(object, random = FALSE, test_log_parms = FALSE, conf.level = 0.6)
    + 
    mean_degparms(
+  object,
+  random = FALSE,
+  test_log_parms = FALSE,
+  conf.level = 0.6,
+  default_log_parms = NA
+)

    Arguments

    object

    An mmkin row object containing several fits of the same model to different datasets

    + +
    random

    Should a list with fixed and random effects be returned?

    + +
    test_log_parms

    If TRUE, log parameters are only considered in the mean calculations if their untransformed counterparts (most likely rate constants) pass the t-test for significant difference from zero.

    + +
    conf.level

    Possibility to adjust the required confidence level for parameter that are tested if requested by 'test_log_parms'.

    + + +
    default_log_parms
    +

    If set to a numeric value, this is used +as a default value for the tested log parameters that failed the +t-test.

    +

    Value

    -

    If random is FALSE (default), a named vector containing mean values + + +

    If random is FALSE (default), a named vector containing mean values of the fitted degradation model parameters. If random is TRUE, a list with fixed and random effects, in the format required by the start argument of nlme for the case of a single grouping variable ds.

    @@ -121,7 +145,7 @@ nlme for the case of a single grouping variable ds.

    -

    Site built with pkgdown 2.0.2.

    +

    Site built with pkgdown 2.0.6.

    diff --git a/docs/dev/reference/mhmkin.html b/docs/dev/reference/mhmkin.html new file mode 100644 index 00000000..e77aace3 --- /dev/null +++ b/docs/dev/reference/mhmkin.html @@ -0,0 +1,213 @@ + +Fit nonlinear mixed-effects models built from one or more kinetic +degradation models and one or more error models — mhmkin • mkin + + +
    +
    + + + +
    +
    +
    +

    Fit nonlinear mixed-effects models built from one or more kinetic +degradation models and one or more error models

    + Source: R/mhmkin.R + +
    + +
    +

    The name of the methods expresses that (multiple) hierarchichal +(also known as multilevel) multicompartment kinetic models are +fitted. Our kinetic models are nonlinear, so we can use various nonlinear +mixed-effects model fitting functions.

    +
    + +
    + 
    mhmkin(objects, backend = "saemix", algorithm = "saem", ...)
+
+# S3 method for list
+mhmkin(
+  objects,
+  backend = "saemix",
+  ...,
+  cores = if (Sys.info()["sysname"] == "Windows") 1 else parallel::detectCores(),
+  cluster = NULL
+)
+
+# S3 method for mhmkin
+[(x, i, j, ..., drop = FALSE)
+
+# S3 method for mhmkin
+print(x, ...)
    +
    + +
    +

    Arguments

    +
    objects
    +

    A list of mmkin objects containing fits of the same +degradation models to the same data, but using different error models.

    + + +
    backend
    +

    The backend to be used for fitting. Currently, only saemix is +supported

    + + +
    algorithm
    +

    The algorithm to be used for fitting (currently not used)

    + + +
    ...
    +

    Further arguments that will be passed to the nonlinear mixed-effects +model fitting function.

    + + +
    cores
    +

    The number of cores to be used for multicore processing. This +is only used when the cluster argument is NULL. On Windows +machines, cores > 1 is not supported, you need to use the cluster +argument to use multiple logical processors. Per default, all cores detected +by parallel::detectCores() are used, except on Windows where the default +is 1.

    + + +
    cluster
    +

    A cluster as returned by makeCluster to be used for +parallel execution.

    + + +
    x
    +

    An mhmkin object.

    + + +
    i
    +

    Row index selecting the fits for specific models

    + + +
    j
    +

    Column index selecting the fits to specific datasets

    + + +
    drop
    +

    If FALSE, the method always returns an mhmkin object, otherwise +either a list of fit objects or a single fit object.

    + +
    +
    +

    Value

    + + +

    A two-dimensional array of fit objects and/or try-errors that can +be indexed using the degradation model names for the first index (row index) +and the error model names for the second index (column index), with class +attribute 'mhmkin'.

    + + +

    An object of class mhmkin.

    +
    +
    +

    See also

    +

    [.mhmkin for subsetting mhmkin objects

    +
    +
    +

    Author

    +

    Johannes Ranke

    +
    + +
    + +
    + + +
    +

    Developed by Johannes Ranke.

    +
    + +
    +

    Site built with pkgdown 2.0.6.

    +
    + +
    + + + + + + + + diff --git a/docs/dev/reference/mkinds.html b/docs/dev/reference/mkinds.html index 0d1de46d..b571e3a0 100644 --- a/docs/dev/reference/mkinds.html +++ b/docs/dev/reference/mkinds.html @@ -1,70 +1,15 @@ - - - - - - - -A dataset class for mkin — mkinds • mkin - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -A dataset class for mkin — mkinds • mkin - - - - - - - - - - - - - + + -
    -
    - -
    - -
    +

    A dataset class for mkin

    - Source: R/mkinds.R + Source: R/mkinds.R
    @@ -153,118 +92,148 @@ such as the on contained in the data field of mkinds objects. Some datasets provided by this package come as mkinds objects nevertheless.

    - 
    # S3 method for mkinds
-print(x, data = FALSE, ...)
    +
    + 
    # S3 method for mkinds
+print(x, data = FALSE, ...)
    +
    + +
    +

    Arguments

    +
    x
    +

    An mkinds object.

    + + +
    data
    +

    Should the data be printed?

    + + +
    ...
    +

    Not used.

    -

    Arguments

    - - - - - - - - - - - - - - -
    x

    An mkinds object.

    data

    Should the data be printed?

    ...

    Not used.

    +
    +
    +

    Public fields

    +

    title
    +

    A full title for the dataset

    -

    Public fields

    -

    -
    title

    A full title for the dataset

    +
    sampling_times
    +

    The sampling times

    -
    sampling_times

    The sampling times

    -
    time_unit

    The time unit

    +
    time_unit
    +

    The time unit

    -
    observed

    Names of the observed variables

    -
    unit

    The unit of the observations

    +
    observed
    +

    Names of the observed variables

    -
    replicates

    The maximum number of replicates per sampling time

    -
    data

    A data frame with at least the columns name, time +

    unit
    +

    The unit of the observations

    + + +
    replicates
    +

    The maximum number of replicates per sampling time

    + + +
    data
    +

    A data frame with at least the columns name, time and value in order to be compatible with mkinfit

    -

    -

    Methods

    +

    +
    +
    +

    Methods

    -

    Public methods

    +
    +

    Public methods

    + +

    +

    Method new()

    +

    Create a new mkinds object

    +

    Usage

    +

    mkinds$new(title = "", data, time_unit = NA, unit = NA)

    +
    - -

    -

    Method new()

    -

    Create a new mkinds object

    Usage

    -

    mkinds$new(title = "", data, time_unit = NA, unit = NA)

    +
    +

    Arguments

    +

    title
    +

    The dataset title

    -

    Arguments

    -

    -
    title

    The dataset title

    -
    data

    The data

    +
    data
    +

    The data

    -
    time_unit

    The time unit

    -
    unit

    The unit of the observations

    +
    time_unit
    +

    The time unit

    -

    -

    -

    Method clone()

    -

    The objects of this class are cloneable with this method.

    Usage

    -

    mkinds$clone(deep = FALSE)

    -

    Arguments

    -

    -
    deep

    Whether to make a deep clone.

    +
    unit
    +

    The unit of the observations

    -

    +

    +
    + +

    +

    Method clone()

    +

    The objects of this class are cloneable with this method.

    +

    Usage

    +

    mkinds$clone(deep = FALSE)

    +
    + +
    +

    Arguments

    +

    deep
    +

    Whether to make a deep clone.

    -

    Examples

    - 
    
-mds <- mkinds$new("FOCUS A", FOCUS_2006_A)
-print(mds)
-
    #> <mkinds> with $title:  FOCUS A 
-#> Observed compounds $observed:  parent 
-#> Sampling times $sampling_times:
-#> 0, 3, 7, 14, 30, 62, 90, 118 
-#> With a maximum of  1  replicates
    
-
    +

    +
    + +
    + +
    + +
    +

    Examples

    + 
    
+mds <- mkinds$new("FOCUS A", FOCUS_2006_A)
+print(mds)
+#> <mkinds> with $title:  FOCUS A 
+#> Observed compounds $observed:  parent 
+#> Sampling times $sampling_times:
+#> 0, 3, 7, 14, 30, 62, 90, 118 
+#> With a maximum of  1  replicates
+
+
    +
    +
    -
    -
    -

    Developed by Johannes Ranke.

    +
    +

    Developed by Johannes Ranke.

    -

    Site built with pkgdown 1.6.1.

    +

    Site built with pkgdown 2.0.6.

    -
    -
    +
    - - + + diff --git a/docs/dev/reference/mkindsg.html b/docs/dev/reference/mkindsg.html index 67c6e5df..d19a7a1d 100644 --- a/docs/dev/reference/mkindsg.html +++ b/docs/dev/reference/mkindsg.html @@ -1,70 +1,15 @@ - - - - - - - -A class for dataset groups for mkin — mkindsg • mkin - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -A class for dataset groups for mkin — mkindsg • mkin - - - - - - - - - - - - - + + -
    -
    - -
    - -
    +

    A class for dataset groups for mkin

    - Source: R/mkinds.R + Source: R/mkinds.R
    @@ -153,308 +92,336 @@ so that combined evaluations are desirable.

    dataset if no data are supplied.

    - 
    # S3 method for mkindsg
-print(x, data = FALSE, verbose = data, ...)
    - -

    Arguments

    - - - - - - - - - - - - - - - - - - -
    x

    An mkindsg object.

    data

    Should the mkinds objects be printed with their data?

    verbose

    Should the mkinds objects be printed?

    ...

    Not used.

    - -

    Public fields

    - -

    -
    title

    A title for the dataset group

    - -
    ds

    A list of mkinds objects

    - -
    observed_n

    Occurrence counts of compounds in datasets

    - -
    f_time_norm

    Time normalisation factors

    - -
    meta

    A data frame with a row for each dataset, +

    + 
    # S3 method for mkindsg
+print(x, data = FALSE, verbose = data, ...)
    +
    + +
    +

    Arguments

    +
    x
    +

    An mkindsg object.

    + + +
    data
    +

    Should the mkinds objects be printed with their data?

    + + +
    verbose
    +

    Should the mkinds objects be printed?

    + + +
    ...
    +

    Not used.

    + +
    +
    +

    Public fields

    +

    title
    +

    A title for the dataset group

    + + +
    ds
    +

    A list of mkinds objects

    + + +
    observed_n
    +

    Occurrence counts of compounds in datasets

    + + +
    f_time_norm
    +

    Time normalisation factors

    + + +
    meta
    +

    A data frame with a row for each dataset, containing additional information in the form of categorical data (factors) or numerical data (e.g. temperature, moisture, or covariates like soil pH).

    -

    -

    Methods

    +

    +
    +
    +

    Methods

    -

    Public methods

    - - -

    -

    Method new()

    -

    Create a new mkindsg object

    Usage

    -

    mkindsg$new(title = "", ds, f_time_norm = rep(1, length(ds)), meta)

    - -

    Arguments

    -

    -
    title

    The title

    - -
    ds

    A list of mkinds objects

    - -
    f_time_norm

    Time normalisation factors

    - -
    meta

    The meta data

    - -

    -

    -

    Method clone()

    -

    The objects of this class are cloneable with this method.

    Usage

    -

    mkindsg$clone(deep = FALSE)

    - -

    Arguments

    -

    -
    deep

    Whether to make a deep clone.

    - -

    - - - -

    Examples

    - 
    
-mdsg <- mkindsg$new("Experimental X", experimental_data_for_UBA_2019[6:10])
-print(mdsg)
-
    #> <mkindsg> holding 5 mkinds objects
-#> Title $title:  Experimental X 
-#> Occurrence of observed compounds $observed_n:
-#> parent     A1 
-#>      5      5 
    print(mdsg, verbose = TRUE)
-
    #> <mkindsg> holding 5 mkinds objects
-#> Title $title:  Experimental X 
-#> Occurrence of observed compounds $observed_n:
-#> parent     A1 
-#>      5      5 
-#> 
-#> Datasets $ds:
-#> <mkinds> with $title:  Soil 6 
-#> Observed compounds $observed:  parent, A1 
-#> Sampling times $sampling_times:
-#> 0, 3, 6, 10, 20, 34, 55, 90, 112, 132 
-#> With a maximum of  2  replicates
-#> Time unit:  days 
-#> Observation unit:  \%AR 
-#> 
-#> <mkinds> with $title:  Soil 7 
-#> Observed compounds $observed:  parent, A1 
-#> Sampling times $sampling_times:
-#> 0, 3, 7, 14, 30, 60, 90, 120, 180 
-#> With a maximum of  2  replicates
-#> Time unit:  days 
-#> Observation unit:  \%AR 
-#> 
-#> <mkinds> with $title:  Soil 8 
-#> Observed compounds $observed:  parent, A1 
-#> Sampling times $sampling_times:
-#> 0, 1, 3, 8, 14, 27, 48, 70 
-#> With a maximum of  2  replicates
-#> Time unit:  days 
-#> Observation unit:  \%AR 
-#> 
-#> <mkinds> with $title:  Soil 9 
-#> Observed compounds $observed:  parent, A1 
-#> Sampling times $sampling_times:
-#> 0, 1, 3, 8, 14, 27, 48, 70, 91, 120 
-#> With a maximum of  2  replicates
-#> Time unit:  days 
-#> Observation unit:  \%AR 
-#> 
-#> <mkinds> with $title:  Soil 10 
-#> Observed compounds $observed:  parent, A1 
-#> Sampling times $sampling_times:
-#> 0, 8, 14, 21, 41, 63, 91, 120 
-#> With a maximum of  2  replicates
-#> Time unit:  days 
-#> Observation unit:  \%AR 
    print(mdsg, verbose = TRUE, data = TRUE)
-
    #> <mkindsg> holding 5 mkinds objects
-#> Title $title:  Experimental X 
-#> Occurrence of observed compounds $observed_n:
-#> parent     A1 
-#>      5      5 
-#> 
-#> Datasets $ds:
-#> <mkinds> with $title:  Soil 6 
-#> Observed compounds $observed:  parent, A1 
-#> Sampling times $sampling_times:
-#> 0, 3, 6, 10, 20, 34, 55, 90, 112, 132 
-#> With a maximum of  2  replicates
-#> Time unit:  days 
-#> Observation unit:  \%AR 
-#>    time parent   A1
-#> 1     0   97.2   NA
-#> 2     0   96.4   NA
-#> 3     3   71.1  4.3
-#> 4     3   69.2  4.6
-#> 5     6   58.1  7.0
-#> 6     6   56.6  7.2
-#> 7    10   44.4  8.2
-#> 8    10   43.4  8.0
-#> 9    20   33.3 11.0
-#> 10   20   29.2 13.7
-#> 11   34   17.6 11.5
-#> 12   34   18.0 12.7
-#> 13   55   10.5 14.9
-#> 14   55    9.3 14.5
-#> 15   90    4.5 12.1
-#> 16   90    4.7 12.3
-#> 17  112    3.0  9.9
-#> 18  112    3.4 10.2
-#> 19  132    2.3  8.8
-#> 20  132    2.7  7.8
-#> 
-#> <mkinds> with $title:  Soil 7 
-#> Observed compounds $observed:  parent, A1 
-#> Sampling times $sampling_times:
-#> 0, 3, 7, 14, 30, 60, 90, 120, 180 
-#> With a maximum of  2  replicates
-#> Time unit:  days 
-#> Observation unit:  \%AR 
-#>    time parent   A1
-#> 1     0   93.6   NA
-#> 2     0   92.3   NA
-#> 3     3   87.0  3.9
-#> 4     3   82.2  3.1
-#> 5     7   74.0  6.9
-#> 6     7   73.9  6.6
-#> 7    14   64.2 10.4
-#> 8    14   69.5  8.3
-#> 9    30   54.0 14.4
-#> 10   30   54.6 13.7
-#> 11   60   41.1 22.1
-#> 12   60   38.4 22.3
-#> 13   90   32.5 27.5
-#> 14   90   35.5 25.4
-#> 15  120   28.1 28.0
-#> 16  120   29.0 26.6
-#> 17  180   26.5 25.8
-#> 18  180   27.6 25.3
-#> 
-#> <mkinds> with $title:  Soil 8 
-#> Observed compounds $observed:  parent, A1 
-#> Sampling times $sampling_times:
-#> 0, 1, 3, 8, 14, 27, 48, 70 
-#> With a maximum of  2  replicates
-#> Time unit:  days 
-#> Observation unit:  \%AR 
-#>    time parent   A1
-#> 1     0   91.9   NA
-#> 2     0   90.8   NA
-#> 3     1   64.9  9.6
-#> 4     1   66.2  7.7
-#> 5     3   43.5 15.0
-#> 6     3   44.1 15.1
-#> 7     8   18.3 21.2
-#> 8     8   18.1 21.1
-#> 9    14   10.2 19.7
-#> 10   14   10.8 18.9
-#> 11   27    4.9 17.5
-#> 12   27    3.3 15.9
-#> 13   48    1.6  9.5
-#> 14   48    1.5  9.8
-#> 15   70    1.1  6.2
-#> 16   70    0.9  6.1
-#> 
-#> <mkinds> with $title:  Soil 9 
-#> Observed compounds $observed:  parent, A1 
-#> Sampling times $sampling_times:
-#> 0, 1, 3, 8, 14, 27, 48, 70, 91, 120 
-#> With a maximum of  2  replicates
-#> Time unit:  days 
-#> Observation unit:  \%AR 
-#>    time parent   A1
-#> 1     0   99.8   NA
-#> 2     0   98.3   NA
-#> 3     1   77.1  4.2
-#> 4     1   77.2  3.9
-#> 5     3   59.0  7.4
-#> 6     3   58.1  7.9
-#> 7     8   27.4 14.5
-#> 8     8   29.2 13.7
-#> 9    14   19.1 14.2
-#> 10   14   29.6 12.2
-#> 11   27   10.1 13.7
-#> 12   27   18.2 13.2
-#> 13   48    4.5 13.6
-#> 14   48    9.1 15.4
-#> 15   70    2.3 10.4
-#> 16   70    2.9 11.6
-#> 17   91    2.0 10.0
-#> 18   91    1.8  9.5
-#> 19  120    2.0  9.1
-#> 20  120    2.2  9.0
-#> 
-#> <mkinds> with $title:  Soil 10 
-#> Observed compounds $observed:  parent, A1 
-#> Sampling times $sampling_times:
-#> 0, 8, 14, 21, 41, 63, 91, 120 
-#> With a maximum of  2  replicates
-#> Time unit:  days 
-#> Observation unit:  \%AR 
-#>    time parent   A1
-#> 1     0   96.1   NA
-#> 2     0   94.3   NA
-#> 3     8   73.9  3.3
-#> 4     8   73.9  3.4
-#> 5    14   69.4  3.9
-#> 6    14   73.1  2.9
-#> 7    21   65.6  6.4
-#> 8    21   65.3  7.2
-#> 9    41   55.9  9.1
-#> 10   41   54.4  8.5
-#> 11   63   47.0 11.7
-#> 12   63   49.3 12.0
-#> 13   91   44.7 13.3
-#> 14   91   46.7 13.2
-#> 15  120   42.1 14.3
-#> 16  120   41.3 12.1
    
-
    +
    +

    Public methods

    + +

    +

    Method new()

    +

    Create a new mkindsg object

    +

    Usage

    +

    mkindsg$new(title = "", ds, f_time_norm = rep(1, length(ds)), meta)

    +
    + +
    +

    Arguments

    +

    title
    +

    The title

    + + +
    ds
    +

    A list of mkinds objects

    + + +
    f_time_norm
    +

    Time normalisation factors

    + + +
    meta
    +

    The meta data

    + + +

    +
    + +

    +

    Method clone()

    +

    The objects of this class are cloneable with this method.

    +

    Usage

    +

    mkindsg$clone(deep = FALSE)

    +
    + +
    +

    Arguments

    +

    deep
    +

    Whether to make a deep clone.

    + + +

    +
    + +
    + +
    + +
    +

    Examples

    + 
    
+mdsg <- mkindsg$new("Experimental X", experimental_data_for_UBA_2019[6:10])
+print(mdsg)
+#> <mkindsg> holding 5 mkinds objects
+#> Title $title:  Experimental X 
+#> Occurrence of observed compounds $observed_n:
+#> parent     A1 
+#>      5      5 
+print(mdsg, verbose = TRUE)
+#> <mkindsg> holding 5 mkinds objects
+#> Title $title:  Experimental X 
+#> Occurrence of observed compounds $observed_n:
+#> parent     A1 
+#>      5      5 
+#> 
+#> Datasets $ds:
+#> <mkinds> with $title:  Soil 6 
+#> Observed compounds $observed:  parent, A1 
+#> Sampling times $sampling_times:
+#> 0, 3, 6, 10, 20, 34, 55, 90, 112, 132 
+#> With a maximum of  2  replicates
+#> Time unit:  days 
+#> Observation unit:  \%AR 
+#> 
+#> <mkinds> with $title:  Soil 7 
+#> Observed compounds $observed:  parent, A1 
+#> Sampling times $sampling_times:
+#> 0, 3, 7, 14, 30, 60, 90, 120, 180 
+#> With a maximum of  2  replicates
+#> Time unit:  days 
+#> Observation unit:  \%AR 
+#> 
+#> <mkinds> with $title:  Soil 8 
+#> Observed compounds $observed:  parent, A1 
+#> Sampling times $sampling_times:
+#> 0, 1, 3, 8, 14, 27, 48, 70 
+#> With a maximum of  2  replicates
+#> Time unit:  days 
+#> Observation unit:  \%AR 
+#> 
+#> <mkinds> with $title:  Soil 9 
+#> Observed compounds $observed:  parent, A1 
+#> Sampling times $sampling_times:
+#> 0, 1, 3, 8, 14, 27, 48, 70, 91, 120 
+#> With a maximum of  2  replicates
+#> Time unit:  days 
+#> Observation unit:  \%AR 
+#> 
+#> <mkinds> with $title:  Soil 10 
+#> Observed compounds $observed:  parent, A1 
+#> Sampling times $sampling_times:
+#> 0, 8, 14, 21, 41, 63, 91, 120 
+#> With a maximum of  2  replicates
+#> Time unit:  days 
+#> Observation unit:  \%AR 
+print(mdsg, verbose = TRUE, data = TRUE)
+#> <mkindsg> holding 5 mkinds objects
+#> Title $title:  Experimental X 
+#> Occurrence of observed compounds $observed_n:
+#> parent     A1 
+#>      5      5 
+#> 
+#> Datasets $ds:
+#> <mkinds> with $title:  Soil 6 
+#> Observed compounds $observed:  parent, A1 
+#> Sampling times $sampling_times:
+#> 0, 3, 6, 10, 20, 34, 55, 90, 112, 132 
+#> With a maximum of  2  replicates
+#> Time unit:  days 
+#> Observation unit:  \%AR 
+#>    time parent   A1
+#> 1     0   97.2   NA
+#> 2     0   96.4   NA
+#> 3     3   71.1  4.3
+#> 4     3   69.2  4.6
+#> 5     6   58.1  7.0
+#> 6     6   56.6  7.2
+#> 7    10   44.4  8.2
+#> 8    10   43.4  8.0
+#> 9    20   33.3 11.0
+#> 10   20   29.2 13.7
+#> 11   34   17.6 11.5
+#> 12   34   18.0 12.7
+#> 13   55   10.5 14.9
+#> 14   55    9.3 14.5
+#> 15   90    4.5 12.1
+#> 16   90    4.7 12.3
+#> 17  112    3.0  9.9
+#> 18  112    3.4 10.2
+#> 19  132    2.3  8.8
+#> 20  132    2.7  7.8
+#> 
+#> <mkinds> with $title:  Soil 7 
+#> Observed compounds $observed:  parent, A1 
+#> Sampling times $sampling_times:
+#> 0, 3, 7, 14, 30, 60, 90, 120, 180 
+#> With a maximum of  2  replicates
+#> Time unit:  days 
+#> Observation unit:  \%AR 
+#>    time parent   A1
+#> 1     0   93.6   NA
+#> 2     0   92.3   NA
+#> 3     3   87.0  3.9
+#> 4     3   82.2  3.1
+#> 5     7   74.0  6.9
+#> 6     7   73.9  6.6
+#> 7    14   64.2 10.4
+#> 8    14   69.5  8.3
+#> 9    30   54.0 14.4
+#> 10   30   54.6 13.7
+#> 11   60   41.1 22.1
+#> 12   60   38.4 22.3
+#> 13   90   32.5 27.5
+#> 14   90   35.5 25.4
+#> 15  120   28.1 28.0
+#> 16  120   29.0 26.6
+#> 17  180   26.5 25.8
+#> 18  180   27.6 25.3
+#> 
+#> <mkinds> with $title:  Soil 8 
+#> Observed compounds $observed:  parent, A1 
+#> Sampling times $sampling_times:
+#> 0, 1, 3, 8, 14, 27, 48, 70 
+#> With a maximum of  2  replicates
+#> Time unit:  days 
+#> Observation unit:  \%AR 
+#>    time parent   A1
+#> 1     0   91.9   NA
+#> 2     0   90.8   NA
+#> 3     1   64.9  9.6
+#> 4     1   66.2  7.7
+#> 5     3   43.5 15.0
+#> 6     3   44.1 15.1
+#> 7     8   18.3 21.2
+#> 8     8   18.1 21.1
+#> 9    14   10.2 19.7
+#> 10   14   10.8 18.9
+#> 11   27    4.9 17.5
+#> 12   27    3.3 15.9
+#> 13   48    1.6  9.5
+#> 14   48    1.5  9.8
+#> 15   70    1.1  6.2
+#> 16   70    0.9  6.1
+#> 
+#> <mkinds> with $title:  Soil 9 
+#> Observed compounds $observed:  parent, A1 
+#> Sampling times $sampling_times:
+#> 0, 1, 3, 8, 14, 27, 48, 70, 91, 120 
+#> With a maximum of  2  replicates
+#> Time unit:  days 
+#> Observation unit:  \%AR 
+#>    time parent   A1
+#> 1     0   99.8   NA
+#> 2     0   98.3   NA
+#> 3     1   77.1  4.2
+#> 4     1   77.2  3.9
+#> 5     3   59.0  7.4
+#> 6     3   58.1  7.9
+#> 7     8   27.4 14.5
+#> 8     8   29.2 13.7
+#> 9    14   19.1 14.2
+#> 10   14   29.6 12.2
+#> 11   27   10.1 13.7
+#> 12   27   18.2 13.2
+#> 13   48    4.5 13.6
+#> 14   48    9.1 15.4
+#> 15   70    2.3 10.4
+#> 16   70    2.9 11.6
+#> 17   91    2.0 10.0
+#> 18   91    1.8  9.5
+#> 19  120    2.0  9.1
+#> 20  120    2.2  9.0
+#> 
+#> <mkinds> with $title:  Soil 10 
+#> Observed compounds $observed:  parent, A1 
+#> Sampling times $sampling_times:
+#> 0, 8, 14, 21, 41, 63, 91, 120 
+#> With a maximum of  2  replicates
+#> Time unit:  days 
+#> Observation unit:  \%AR 
+#>    time parent   A1
+#> 1     0   96.1   NA
+#> 2     0   94.3   NA
+#> 3     8   73.9  3.3
+#> 4     8   73.9  3.4
+#> 5    14   69.4  3.9
+#> 6    14   73.1  2.9
+#> 7    21   65.6  6.4
+#> 8    21   65.3  7.2
+#> 9    41   55.9  9.1
+#> 10   41   54.4  8.5
+#> 11   63   47.0 11.7
+#> 12   63   49.3 12.0
+#> 13   91   44.7 13.3
+#> 14   91   46.7 13.2
+#> 15  120   42.1 14.3
+#> 16  120   41.3 12.1
+
+
    +
    +
    -
    -
    -

    Developed by Johannes Ranke.

    +
    +

    Developed by Johannes Ranke.

    -

    Site built with pkgdown 1.6.1.

    +

    Site built with pkgdown 2.0.6.

    -
    -
    +
    - - + + diff --git a/docs/dev/reference/mkinfit-1.png b/docs/dev/reference/mkinfit-1.png index de2a90a9..7c51deb6 100644 Binary files a/docs/dev/reference/mkinfit-1.png and b/docs/dev/reference/mkinfit-1.png differ diff --git a/docs/dev/reference/mkinfit.html b/docs/dev/reference/mkinfit.html index 5910038f..a00a0458 100644 --- a/docs/dev/reference/mkinfit.html +++ b/docs/dev/reference/mkinfit.html @@ -1,46 +1,5 @@ - - - - - - - -Fit a kinetic model to data with one or more state variables — mkinfit • mkin - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -Fit a kinetic model to data with one or more state variables — mkinfit • mkin - - - - - - - - - + + - - - - -
    -
    - -
    - -
    +

    Fit a kinetic model to data with one or more state variables

    - Source: R/mkinfit.R + Source: R/mkinfit.R

    This function maximises the likelihood of the observed data using the Port -algorithm stats::nlminb(), and the specified initial or fixed +algorithm stats::nlminb(), and the specified initial or fixed parameters and starting values. In each step of the optimisation, the -kinetic model is solved using the function mkinpredict(), except +kinetic model is solved using the function mkinpredict(), except if an analytical solution is implemented, in which case the model is solved -using the degradation function in the mkinmod object. The +using the degradation function in the mkinmod object. The parameters of the selected error model are fitted simultaneously with the degradation model parameters, as both of them are arguments of the likelihood function.

    - 
    mkinfit(
-  mkinmod,
-  observed,
-  parms.ini = "auto",
-  state.ini = "auto",
-  err.ini = "auto",
-  fixed_parms = NULL,
-  fixed_initials = names(mkinmod$diffs)[-1],
-  from_max_mean = FALSE,
-  solution_type = c("auto", "analytical", "eigen", "deSolve"),
-  method.ode = "lsoda",
-  use_compiled = "auto",
-  control = list(eval.max = 300, iter.max = 200),
-  transform_rates = TRUE,
-  transform_fractions = TRUE,
-  quiet = FALSE,
-  atol = 1e-08,
-  rtol = 1e-10,
-  error_model = c("const", "obs", "tc"),
-  error_model_algorithm = c("auto", "d_3", "direct", "twostep", "threestep",
-    "fourstep", "IRLS", "OLS"),
-  reweight.tol = 1e-08,
-  reweight.max.iter = 10,
-  trace_parms = FALSE,
-  test_residuals = FALSE,
-  ...
-)
    - -

    Arguments

    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    mkinmod

    A list of class mkinmod, containing the kinetic +

    + 
    mkinfit(
+  mkinmod,
+  observed,
+  parms.ini = "auto",
+  state.ini = "auto",
+  err.ini = "auto",
+  fixed_parms = NULL,
+  fixed_initials = names(mkinmod$diffs)[-1],
+  from_max_mean = FALSE,
+  solution_type = c("auto", "analytical", "eigen", "deSolve"),
+  method.ode = "lsoda",
+  use_compiled = "auto",
+  control = list(eval.max = 300, iter.max = 200),
+  transform_rates = TRUE,
+  transform_fractions = TRUE,
+  quiet = FALSE,
+  atol = 1e-08,
+  rtol = 1e-10,
+  error_model = c("const", "obs", "tc"),
+  error_model_algorithm = c("auto", "d_3", "direct", "twostep", "threestep", "fourstep",
+    "IRLS", "OLS"),
+  reweight.tol = 1e-08,
+  reweight.max.iter = 10,
+  trace_parms = FALSE,
+  test_residuals = FALSE,
+  ...
+)
    +
    + +
    +

    Arguments

    +
    mkinmod
    +

    A list of class mkinmod, containing the kinetic model to be fitted to the data, or one of the shorthand names ("SFO", "FOMC", "DFOP", "HS", "SFORB", "IORE"). If a shorthand name is given, a parent only degradation model is generated for the variable with the -highest value in observed.

    observed

    A dataframe with the observed data. The first column called +highest value in observed.

    + + +
    observed
    +

    A dataframe with the observed data. The first column called "name" must contain the name of the observed variable for each data point. The second column must contain the times of observation, named "time". The third column must be named "value" and contain the observed values. @@ -212,11 +151,11 @@ Zero values in the "value" column will be removed, with a warning, in order to avoid problems with fitting the two-component error model. This is not expected to be a problem, because in general, values of zero are not observed in degradation data, because there is a lower limit of -detection.

    parms.ini

    A named vector of initial values for the parameters, +detection.

    + + +
    parms.ini
    +

    A named vector of initial values for the parameters, including parameters to be optimised and potentially also fixed parameters as indicated by fixed_parms. If set to "auto", initial values for rate constants are set to default values. Using parameter names that are @@ -225,113 +164,113 @@ not in the model gives an error.

    needs. You can use the parameter lists "bparms.ode" from a previously fitted model, which contains the differential equation parameters from this model. This works nicely if the models are nested. An example is -given below.

    state.ini

    A named vector of initial values for the state variables of +given below.

    + + +
    state.ini
    +

    A named vector of initial values for the state variables of the model. In case the observed variables are represented by more than one model variable, the names will differ from the names of the observed -variables (see map component of mkinmod). The default +variables (see map component of mkinmod). The default is to set the initial value of the first model variable to the mean of the time zero values for the variable with the maximum observed value, and all others to 0. If this variable has no time zero observations, its initial -value is set to 100.

    err.ini

    A named vector of initial values for the error model +value is set to 100.

    + + +
    err.ini
    +

    A named vector of initial values for the error model parameters to be optimised. If set to "auto", initial values are set to default values. Otherwise, inital values for all error model parameters -must be given.

    fixed_parms

    The names of parameters that should not be optimised but +must be given.

    + + +
    fixed_parms
    +

    The names of parameters that should not be optimised but rather kept at the values specified in parms.ini. Alternatively, a named numeric vector of parameters to be fixed, regardless of the values -in parms.ini.

    fixed_initials

    The names of model variables for which the initial +in parms.ini.

    + + +
    fixed_initials
    +

    The names of model variables for which the initial state at time 0 should be excluded from the optimisation. Defaults to all -state variables except for the first one.

    from_max_mean

    If this is set to TRUE, and the model has only one +state variables except for the first one.

    + + +
    from_max_mean
    +

    If this is set to TRUE, and the model has only one observed variable, then data before the time of the maximum observed value (after averaging for each sampling time) are discarded, and this time is subtracted from all remaining time values, so the time of the maximum -observed mean value is the new time zero.

    solution_type

    If set to "eigen", the solution of the system of +observed mean value is the new time zero.

    + + +
    solution_type
    +

    If set to "eigen", the solution of the system of differential equations is based on the spectral decomposition of the coefficient matrix in cases that this is possible. If set to "deSolve", a -numerical ode solver from package deSolve is used. If +numerical ode solver from package deSolve is used. If set to "analytical", an analytical solution of the model is used. This is only implemented for relatively simple degradation models. The default is "auto", which uses "analytical" if possible, otherwise "deSolve" if a compiler is present, and "eigen" if no compiler is present and the model -can be expressed using eigenvalues and eigenvectors.

    method.ode

    The solution method passed via mkinpredict() -to deSolve::ode() in case the solution type is "deSolve". The default -"lsoda" is performant, but sometimes fails to converge.

    use_compiled

    If set to FALSE, no compiled version of the -mkinmod model is used in the calls to mkinpredict() even if a compiled -version is present.

    control

    A list of control arguments passed to stats::nlminb().

    transform_rates

    Boolean specifying if kinetic rate constants should +can be expressed using eigenvalues and eigenvectors.

    + + +
    method.ode
    +

    The solution method passed via mkinpredict() +to deSolve::ode() in case the solution type is "deSolve". The default +"lsoda" is performant, but sometimes fails to converge.

    + + +
    use_compiled
    +

    If set to FALSE, no compiled version of the +mkinmod model is used in the calls to mkinpredict() even if a compiled +version is present.

    + + +
    control
    +

    A list of control arguments passed to stats::nlminb().

    + + +
    transform_rates
    +

    Boolean specifying if kinetic rate constants should be transformed in the model specification used in the fitting for better compliance with the assumption of normal distribution of the estimator. If TRUE, also alpha and beta parameters of the FOMC model are log-transformed, as well as k1 and k2 rate constants for the DFOP and HS models and the break point tb of the HS model. If FALSE, zero is used as -a lower bound for the rates in the optimisation.

    transform_fractions

    Boolean specifying if formation fractions +a lower bound for the rates in the optimisation.

    + + +
    transform_fractions
    +

    Boolean specifying if formation fractions should be transformed in the model specification used in the fitting for better compliance with the assumption of normal distribution of the estimator. The default (TRUE) is to do transformations. If TRUE, the g parameter of the DFOP model is also transformed. Transformations -are described in transform_odeparms.

    quiet

    Suppress printing out the current value of the negative -log-likelihood after each improvement?

    atol

    Absolute error tolerance, passed to deSolve::ode(). Default -is 1e-8, which is lower than the default in the deSolve::lsoda() -function which is used per default.

    rtol

    Absolute error tolerance, passed to deSolve::ode(). Default -is 1e-10, much lower than in deSolve::lsoda().

    error_model

    If the error model is "const", a constant standard +are described in transform_odeparms.

    + + +
    quiet
    +

    Suppress printing out the current value of the negative +log-likelihood after each improvement?

    + + +
    atol
    +

    Absolute error tolerance, passed to deSolve::ode(). Default +is 1e-8, which is lower than the default in the deSolve::lsoda() +function which is used per default.

    + + +
    rtol
    +

    Absolute error tolerance, passed to deSolve::ode(). Default +is 1e-10, much lower than in deSolve::lsoda().

    + + +
    error_model
    +

    If the error model is "const", a constant standard deviation is assumed.

    If the error model is "obs", each observed variable is assumed to have its own variance.

    @@ -340,11 +279,11 @@ error model similar to the one described by Rocke and Lorenzato (1995) is used for setting up the likelihood function. Note that this model deviates from the model by Rocke and Lorenzato, as their model implies that the errors follow a lognormal distribution for large values, not a -normal distribution as assumed by this method.

    error_model_algorithm

    If "auto", the selected algorithm depends on +normal distribution as assumed by this method.

    + + +
    error_model_algorithm
    +

    If "auto", the selected algorithm depends on the error model. If the error model is "const", unweighted nonlinear least squares fitting ("OLS") is selected. If the error model is "obs", or "tc", the "d_3" algorithm is selected.

    @@ -367,372 +306,391 @@ degradation and error model parameters.

    unweighted least squares, and then iterates optimization of the error model parameters and subsequent optimization of the degradation model using those error model parameters, until the error model parameters -converge.

    reweight.tol

    Tolerance for the convergence criterion calculated from -the error model parameters in IRLS fits.

    reweight.max.iter

    Maximum number of iterations in IRLS fits.

    trace_parms

    Should a trace of the parameter values be listed?

    test_residuals

    Should the residuals be tested for normal distribution?

    ...

    Further arguments that will be passed on to -deSolve::ode().

    - -

    Value

    - -

    A list with "mkinfit" in the class attribute.

    -

    Details

    +converge.

    + + +
    reweight.tol
    +

    Tolerance for the convergence criterion calculated from +the error model parameters in IRLS fits.

    + +
    reweight.max.iter
    +

    Maximum number of iterations in IRLS fits.

    + + +
    trace_parms
    +

    Should a trace of the parameter values be listed?

    + + +
    test_residuals
    +

    Should the residuals be tested for normal distribution?

    + + +
    ...
    +

    Further arguments that will be passed on to +deSolve::ode().

    + +
    +
    +

    Value

    + + +

    A list with "mkinfit" in the class attribute.

    +
    +
    +

    Details

    Per default, parameters in the kinetic models are internally transformed in order to better satisfy the assumption of a normal distribution of their estimators.

    -

    Note

    - +
    +
    +

    Note

    When using the "IORE" submodel for metabolites, fitting with "transform_rates = TRUE" (the default) often leads to failures of the numerical ODE solver. In this situation it may help to switch off the internal rate transformation.

    -

    References

    - +
    +
    +

    References

    Rocke DM and Lorenzato S (1995) A two-component model for measurement error in analytical chemistry. Technometrics 37(2), 176-184.

    Ranke J and Meinecke S (2019) Error Models for the Kinetic Evaluation of Chemical Degradation Data. Environments 6(12) 124 -doi: 10.3390/environments6120124 +doi:10.3390/environments6120124 .

    -

    See also

    - -

    summary.mkinfit, plot.mkinfit, parms and lrtest.

    +
    +
    +

    See also

    +

    summary.mkinfit, plot.mkinfit, parms and lrtest.

    Comparisons of models fitted to the same data can be made using -AIC by virtue of the method logLik.mkinfit.

    +AIC by virtue of the method logLik.mkinfit.

    Fitting of several models to several datasets in a single call to -mmkin.

    -

    Author

    - +mmkin.

    +
    +
    +

    Author

    Johannes Ranke

    +
    -

    Examples

    - 
    
-# Use shorthand notation for parent only degradation
-fit <- mkinfit("FOMC", FOCUS_2006_C, quiet = TRUE)
-summary(fit)
-
    #> mkin version used for fitting:    1.0.3.9000 
-#> R version used for fitting:       4.0.3 
-#> Date of fit:     Mon Feb 15 17:09:39 2021 
-#> Date of summary: Mon Feb 15 17:09:39 2021 
-#> 
-#> Equations:
-#> d_parent/dt = - (alpha/beta) * 1/((time/beta) + 1) * parent
-#> 
-#> Model predictions using solution type analytical 
-#> 
-#> Fitted using 222 model solutions performed in 0.045 s
-#> 
-#> Error model: Constant variance 
-#> 
-#> Error model algorithm: OLS 
-#> 
-#> Starting values for parameters to be optimised:
-#>          value   type
-#> parent_0  85.1  state
-#> alpha      1.0 deparm
-#> beta      10.0 deparm
-#> 
-#> Starting values for the transformed parameters actually optimised:
-#>               value lower upper
-#> parent_0  85.100000  -Inf   Inf
-#> log_alpha  0.000000  -Inf   Inf
-#> log_beta   2.302585  -Inf   Inf
-#> 
-#> Fixed parameter values:
-#> None
-#> 
-#> Results:
-#> 
-#>        AIC      BIC    logLik
-#>   44.68652 45.47542 -18.34326
-#> 
-#> Optimised, transformed parameters with symmetric confidence intervals:
-#>           Estimate Std. Error    Lower   Upper
-#> parent_0  85.87000     1.8070 81.23000 90.5200
-#> log_alpha  0.05192     0.1353 -0.29580  0.3996
-#> log_beta   0.65100     0.2287  0.06315  1.2390
-#> sigma      1.85700     0.4378  0.73200  2.9830
-#> 
-#> Parameter correlation:
-#>             parent_0  log_alpha   log_beta     sigma
-#> parent_0   1.000e+00 -1.565e-01 -3.142e-01 4.772e-08
-#> log_alpha -1.565e-01  1.000e+00  9.564e-01 1.005e-07
-#> log_beta  -3.142e-01  9.564e-01  1.000e+00 8.541e-08
-#> sigma      4.772e-08  1.005e-07  8.541e-08 1.000e+00
-#> 
-#> Backtransformed parameters:
-#> Confidence intervals for internally transformed parameters are asymmetric.
-#> t-test (unrealistically) based on the assumption of normal distribution
-#> for estimators of untransformed parameters.
-#>          Estimate t value    Pr(>t)   Lower  Upper
-#> parent_0   85.870  47.530 3.893e-08 81.2300 90.520
-#> alpha       1.053   7.393 3.562e-04  0.7439  1.491
-#> beta        1.917   4.373 3.601e-03  1.0650  3.451
-#> sigma       1.857   4.243 4.074e-03  0.7320  2.983
-#> 
-#> FOCUS Chi2 error levels in percent:
-#>          err.min n.optim df
-#> All data   6.657       3  6
-#> parent     6.657       3  6
-#> 
-#> Estimated disappearance times:
-#>         DT50  DT90 DT50back
-#> parent 1.785 15.15     4.56
-#> 
-#> Data:
-#>  time variable observed predicted residual
-#>     0   parent     85.1    85.875  -0.7749
-#>     1   parent     57.9    55.191   2.7091
-#>     3   parent     29.9    31.845  -1.9452
-#>     7   parent     14.6    17.012  -2.4124
-#>    14   parent      9.7     9.241   0.4590
-#>    28   parent      6.6     4.754   1.8460
-#>    63   parent      4.0     2.102   1.8977
-#>    91   parent      3.9     1.441   2.4590
-#>   119   parent      0.6     1.092  -0.4919
    
-# One parent compound, one metabolite, both single first order.
-# We remove zero values from FOCUS dataset D in order to avoid warnings
-FOCUS_D <- subset(FOCUS_2006_D, value != 0)
-# Use mkinsub for convenience in model formulation. Pathway to sink included per default.
-SFO_SFO <- mkinmod(
-  parent = mkinsub("SFO", "m1"),
-  m1 = mkinsub("SFO"))
-
    #> Temporary DLL for differentials generated and loaded
    
-# Fit the model quietly to the FOCUS example dataset D using defaults
-fit <- mkinfit(SFO_SFO, FOCUS_D, quiet = TRUE)
-plot_sep(fit)
-
    # As lower parent values appear to have lower variance, we try an alternative error model
-fit.tc <- mkinfit(SFO_SFO, FOCUS_D, quiet = TRUE, error_model = "tc")
-# This avoids the warning, and the likelihood ratio test confirms it is preferable
-lrtest(fit.tc, fit)
-
    #> Likelihood ratio test
-#> 
-#> Model 1: SFO_SFO with error model tc and fixed parameter(s) m1_0
-#> Model 2: SFO_SFO with error model const and fixed parameter(s) m1_0
-#>   #Df  LogLik Df  Chisq Pr(>Chisq)    
-#> 1   6 -64.983                         
-#> 2   5 -97.224 -1 64.483  9.737e-16 ***
-#> ---
-#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
    # We can also allow for different variances of parent and metabolite as error model
-fit.obs <- mkinfit(SFO_SFO, FOCUS_D, quiet = TRUE, error_model = "obs")
-# The two-component error model has significantly higher likelihood
-lrtest(fit.obs, fit.tc)
-
    #> Likelihood ratio test
-#> 
-#> Model 1: SFO_SFO with error model tc and fixed parameter(s) m1_0
-#> Model 2: SFO_SFO with error model obs and fixed parameter(s) m1_0
-#>   #Df  LogLik Df  Chisq Pr(>Chisq)    
-#> 1   6 -64.983                         
-#> 2   6 -96.936  0 63.907  < 2.2e-16 ***
-#> ---
-#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
    parms(fit.tc)
-
    #>       parent_0       k_parent           k_m1 f_parent_to_m1      sigma_low 
-#>   1.007343e+02   1.005562e-01   5.166712e-03   5.083933e-01   3.049883e-03 
-#>       rsd_high 
-#>   7.928118e-02 
    endpoints(fit.tc)
-
    #> $ff
-#>   parent_m1 parent_sink 
-#>   0.5083933   0.4916067 
-#> 
-#> $distimes
-#>             DT50      DT90
-#> parent   6.89313  22.89848
-#> m1     134.15634 445.65772
-#> 
    
-# We can show a quick (only one replication) benchmark for this case, as we
-# have several alternative solution methods for the model. We skip
-# uncompiled deSolve, as it is so slow. More benchmarks are found in the
-# benchmark vignette
-# \dontrun{
-if(require(rbenchmark)) {
-  benchmark(replications = 1, order = "relative", columns = c("test", "relative", "elapsed"),
-    deSolve_compiled = mkinfit(SFO_SFO, FOCUS_D, quiet = TRUE, error_model = "tc",
-      solution_type = "deSolve", use_compiled = TRUE),
-    eigen = mkinfit(SFO_SFO, FOCUS_D, quiet = TRUE, error_model = "tc",
-      solution_type = "eigen"),
-    analytical = mkinfit(SFO_SFO, FOCUS_D, quiet = TRUE, error_model = "tc",
-      solution_type = "analytical"))
-}
-
    #>               test relative elapsed
-#> 3       analytical    1.000   0.563
-#> 1 deSolve_compiled    1.702   0.958
-#> 2            eigen    2.597   1.462
    # }
-
-# Use stepwise fitting, using optimised parameters from parent only fit, FOMC-SFO
-# \dontrun{
-FOMC_SFO <- mkinmod(
-  parent = mkinsub("FOMC", "m1"),
-  m1 = mkinsub("SFO"))
-
    #> Temporary DLL for differentials generated and loaded
    fit.FOMC_SFO <- mkinfit(FOMC_SFO, FOCUS_D, quiet = TRUE)
-# Again, we get a warning and try a more sophisticated error model
-fit.FOMC_SFO.tc <- mkinfit(FOMC_SFO, FOCUS_D, quiet = TRUE, error_model = "tc")
-# This model has a higher likelihood, but not significantly so
-lrtest(fit.tc, fit.FOMC_SFO.tc)
-
    #> Likelihood ratio test
-#> 
-#> Model 1: FOMC_SFO with error model tc and fixed parameter(s) m1_0
-#> Model 2: SFO_SFO with error model tc and fixed parameter(s) m1_0
-#>   #Df  LogLik Df  Chisq Pr(>Chisq)
-#> 1   7 -64.829                     
-#> 2   6 -64.983 -1 0.3075     0.5792
    # Also, the missing standard error for log_beta and the t-tests for alpha
-# and beta indicate overparameterisation
-summary(fit.FOMC_SFO.tc, data = FALSE)
-
    #> Warning: NaNs produced
    #> Warning: NaNs produced
    #> Warning: diag(.) had 0 or NA entries; non-finite result is doubtful
    #> mkin version used for fitting:    1.0.3.9000 
-#> R version used for fitting:       4.0.3 
-#> Date of fit:     Mon Feb 15 17:09:50 2021 
-#> Date of summary: Mon Feb 15 17:09:50 2021 
-#> 
-#> Equations:
-#> d_parent/dt = - (alpha/beta) * 1/((time/beta) + 1) * parent
-#> d_m1/dt = + f_parent_to_m1 * (alpha/beta) * 1/((time/beta) + 1) *
-#>            parent - k_m1 * m1
-#> 
-#> Model predictions using solution type deSolve 
-#> 
-#> Fitted using 3729 model solutions performed in 2.815 s
-#> 
-#> Error model: Two-component variance function 
-#> 
-#> Error model algorithm: d_3 
-#> Direct fitting and three-step fitting yield approximately the same likelihood 
-#> 
-#> Starting values for parameters to be optimised:
-#>                 value   type
-#> parent_0       100.75  state
-#> alpha            1.00 deparm
-#> beta            10.00 deparm
-#> k_m1             0.10 deparm
-#> f_parent_to_m1   0.50 deparm
-#> sigma_low        0.10  error
-#> rsd_high         0.10  error
-#> 
-#> Starting values for the transformed parameters actually optimised:
-#>                      value lower upper
-#> parent_0        100.750000  -Inf   Inf
-#> log_k_m1         -2.302585  -Inf   Inf
-#> f_parent_qlogis   0.000000  -Inf   Inf
-#> log_alpha         0.000000  -Inf   Inf
-#> log_beta          2.302585  -Inf   Inf
-#> sigma_low         0.100000     0   Inf
-#> rsd_high          0.100000     0   Inf
-#> 
-#> Fixed parameter values:
-#>      value  type
-#> m1_0     0 state
-#> 
-#> Results:
-#> 
-#>       AIC      BIC    logLik
-#>   143.658 155.1211 -64.82902
-#> 
-#> Optimised, transformed parameters with symmetric confidence intervals:
-#>                   Estimate Std. Error     Lower      Upper
-#> parent_0        101.600000  2.6400000 96.240000 107.000000
-#> log_k_m1         -5.284000  0.0929100 -5.474000  -5.095000
-#> f_parent_qlogis   0.001426  0.0767000 -0.155000   0.157800
-#> log_alpha         5.522000  0.0077320  5.506000   5.538000
-#> log_beta          7.806000        NaN       NaN        NaN
-#> sigma_low         0.002488  0.0002431  0.001992   0.002984
-#> rsd_high          0.079210  0.0093280  0.060180   0.098230
-#> 
-#> Parameter correlation:
-#>                  parent_0  log_k_m1 f_parent_qlogis log_alpha log_beta
-#> parent_0         1.000000 -0.095226        -0.76678   0.70544      NaN
-#> log_k_m1        -0.095226  1.000000         0.51432  -0.14387      NaN
-#> f_parent_qlogis -0.766780  0.514321         1.00000  -0.61396      NaN
-#> log_alpha        0.705444 -0.143872        -0.61396   1.00000      NaN
-#> log_beta              NaN       NaN             NaN       NaN        1
-#> sigma_low        0.016073  0.001586         0.01548   5.87007      NaN
-#> rsd_high         0.006626 -0.011700        -0.05357   0.04849      NaN
-#>                 sigma_low  rsd_high
-#> parent_0         0.016073  0.006626
-#> log_k_m1         0.001586 -0.011700
-#> f_parent_qlogis  0.015476 -0.053566
-#> log_alpha        5.870075  0.048487
-#> log_beta              NaN       NaN
-#> sigma_low        1.000000 -0.652558
-#> rsd_high        -0.652558  1.000000
-#> 
-#> Backtransformed parameters:
-#> Confidence intervals for internally transformed parameters are asymmetric.
-#> t-test (unrealistically) based on the assumption of normal distribution
-#> for estimators of untransformed parameters.
-#>                 Estimate t value    Pr(>t)     Lower     Upper
-#> parent_0       1.016e+02 32.7800 6.311e-26 9.624e+01 1.070e+02
-#> k_m1           5.072e-03 10.1200 1.216e-11 4.196e-03 6.130e-03
-#> f_parent_to_m1 5.004e-01 20.8300 4.317e-20 4.613e-01 5.394e-01
-#> alpha          2.502e+02  0.5624 2.889e-01 2.463e+02 2.542e+02
-#> beta           2.455e+03  0.5549 2.915e-01        NA        NA
-#> sigma_low      2.488e-03  0.4843 3.158e-01 1.992e-03 2.984e-03
-#> rsd_high       7.921e-02  8.4300 8.001e-10 6.018e-02 9.823e-02
-#> 
-#> FOCUS Chi2 error levels in percent:
-#>          err.min n.optim df
-#> All data   6.781       5 14
-#> parent     7.141       3  6
-#> m1         4.640       2  8
-#> 
-#> Resulting formation fractions:
-#>                 ff
-#> parent_m1   0.5004
-#> parent_sink 0.4996
-#> 
-#> Estimated disappearance times:
-#>           DT50  DT90 DT50back
-#> parent   6.812  22.7    6.834
-#> m1     136.661 454.0       NA
    
-# We can easily use starting parameters from the parent only fit (only for illustration)
-fit.FOMC = mkinfit("FOMC", FOCUS_2006_D, quiet = TRUE, error_model = "tc")
-fit.FOMC_SFO <- mkinfit(FOMC_SFO, FOCUS_D, quiet = TRUE,
-  parms.ini = fit.FOMC$bparms.ode, error_model = "tc")
-# }
-
    +
    +

    Examples

    + 
    
+# Use shorthand notation for parent only degradation
+fit <- mkinfit("FOMC", FOCUS_2006_C, quiet = TRUE)
+summary(fit)
+#> mkin version used for fitting:    1.1.2 
+#> R version used for fitting:       4.2.1 
+#> Date of fit:     Wed Aug 10 15:24:16 2022 
+#> Date of summary: Wed Aug 10 15:24:16 2022 
+#> 
+#> Equations:
+#> d_parent/dt = - (alpha/beta) * 1/((time/beta) + 1) * parent
+#> 
+#> Model predictions using solution type analytical 
+#> 
+#> Fitted using 222 model solutions performed in 0.045 s
+#> 
+#> Error model: Constant variance 
+#> 
+#> Error model algorithm: OLS 
+#> 
+#> Starting values for parameters to be optimised:
+#>          value   type
+#> parent_0  85.1  state
+#> alpha      1.0 deparm
+#> beta      10.0 deparm
+#> 
+#> Starting values for the transformed parameters actually optimised:
+#>               value lower upper
+#> parent_0  85.100000  -Inf   Inf
+#> log_alpha  0.000000  -Inf   Inf
+#> log_beta   2.302585  -Inf   Inf
+#> 
+#> Fixed parameter values:
+#> None
+#> 
+#> Results:
+#> 
+#>        AIC      BIC    logLik
+#>   44.68652 45.47542 -18.34326
+#> 
+#> Optimised, transformed parameters with symmetric confidence intervals:
+#>           Estimate Std. Error    Lower   Upper
+#> parent_0  85.87000     1.8070 81.23000 90.5200
+#> log_alpha  0.05192     0.1353 -0.29580  0.3996
+#> log_beta   0.65100     0.2287  0.06315  1.2390
+#> sigma      1.85700     0.4378  0.73200  2.9830
+#> 
+#> Parameter correlation:
+#>             parent_0  log_alpha   log_beta     sigma
+#> parent_0   1.000e+00 -1.565e-01 -3.142e-01 4.772e-08
+#> log_alpha -1.565e-01  1.000e+00  9.564e-01 1.005e-07
+#> log_beta  -3.142e-01  9.564e-01  1.000e+00 8.541e-08
+#> sigma      4.772e-08  1.005e-07  8.541e-08 1.000e+00
+#> 
+#> Backtransformed parameters:
+#> Confidence intervals for internally transformed parameters are asymmetric.
+#> t-test (unrealistically) based on the assumption of normal distribution
+#> for estimators of untransformed parameters.
+#>          Estimate t value    Pr(>t)   Lower  Upper
+#> parent_0   85.870  47.530 3.893e-08 81.2300 90.520
+#> alpha       1.053   7.393 3.562e-04  0.7439  1.491
+#> beta        1.917   4.373 3.601e-03  1.0650  3.451
+#> sigma       1.857   4.243 4.074e-03  0.7320  2.983
+#> 
+#> FOCUS Chi2 error levels in percent:
+#>          err.min n.optim df
+#> All data   6.657       3  6
+#> parent     6.657       3  6
+#> 
+#> Estimated disappearance times:
+#>         DT50  DT90 DT50back
+#> parent 1.785 15.15     4.56
+#> 
+#> Data:
+#>  time variable observed predicted residual
+#>     0   parent     85.1    85.875  -0.7749
+#>     1   parent     57.9    55.191   2.7091
+#>     3   parent     29.9    31.845  -1.9452
+#>     7   parent     14.6    17.012  -2.4124
+#>    14   parent      9.7     9.241   0.4590
+#>    28   parent      6.6     4.754   1.8460
+#>    63   parent      4.0     2.102   1.8977
+#>    91   parent      3.9     1.441   2.4590
+#>   119   parent      0.6     1.092  -0.4919
+
+# One parent compound, one metabolite, both single first order.
+# We remove zero values from FOCUS dataset D in order to avoid warnings
+FOCUS_D <- subset(FOCUS_2006_D, value != 0)
+# Use mkinsub for convenience in model formulation. Pathway to sink included per default.
+SFO_SFO <- mkinmod(
+  parent = mkinsub("SFO", "m1"),
+  m1 = mkinsub("SFO"))
+#> Temporary DLL for differentials generated and loaded
+
+# Fit the model quietly to the FOCUS example dataset D using defaults
+fit <- mkinfit(SFO_SFO, FOCUS_D, quiet = TRUE)
+plot_sep(fit)
+
+# As lower parent values appear to have lower variance, we try an alternative error model
+fit.tc <- mkinfit(SFO_SFO, FOCUS_D, quiet = TRUE, error_model = "tc")
+# This avoids the warning, and the likelihood ratio test confirms it is preferable
+lrtest(fit.tc, fit)
+#> Likelihood ratio test
+#> 
+#> Model 1: SFO_SFO with error model tc and fixed parameter(s) m1_0
+#> Model 2: SFO_SFO with error model const and fixed parameter(s) m1_0
+#>   #Df  LogLik Df  Chisq Pr(>Chisq)    
+#> 1   6 -64.983                         
+#> 2   5 -97.224 -1 64.483  9.737e-16 ***
+#> ---
+#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+# We can also allow for different variances of parent and metabolite as error model
+fit.obs <- mkinfit(SFO_SFO, FOCUS_D, quiet = TRUE, error_model = "obs")
+# The two-component error model has significantly higher likelihood
+lrtest(fit.obs, fit.tc)
+#> Likelihood ratio test
+#> 
+#> Model 1: SFO_SFO with error model tc and fixed parameter(s) m1_0
+#> Model 2: SFO_SFO with error model obs and fixed parameter(s) m1_0
+#>   #Df  LogLik Df  Chisq Pr(>Chisq)    
+#> 1   6 -64.983                         
+#> 2   6 -96.936  0 63.907  < 2.2e-16 ***
+#> ---
+#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
+parms(fit.tc)
+#>       parent_0       k_parent           k_m1 f_parent_to_m1      sigma_low 
+#>   1.007343e+02   1.005562e-01   5.166712e-03   5.083933e-01   3.049883e-03 
+#>       rsd_high 
+#>   7.928118e-02 
+endpoints(fit.tc)
+#> $ff
+#>   parent_m1 parent_sink 
+#>   0.5083933   0.4916067 
+#> 
+#> $distimes
+#>             DT50      DT90
+#> parent   6.89313  22.89848
+#> m1     134.15634 445.65772
+#> 
+
+# We can show a quick (only one replication) benchmark for this case, as we
+# have several alternative solution methods for the model. We skip
+# uncompiled deSolve, as it is so slow. More benchmarks are found in the
+# benchmark vignette
+# \dontrun{
+if(require(rbenchmark)) {
+  benchmark(replications = 1, order = "relative", columns = c("test", "relative", "elapsed"),
+    deSolve_compiled = mkinfit(SFO_SFO, FOCUS_D, quiet = TRUE, error_model = "tc",
+      solution_type = "deSolve", use_compiled = TRUE),
+    eigen = mkinfit(SFO_SFO, FOCUS_D, quiet = TRUE, error_model = "tc",
+      solution_type = "eigen"),
+    analytical = mkinfit(SFO_SFO, FOCUS_D, quiet = TRUE, error_model = "tc",
+      solution_type = "analytical"))
+}
+#>               test relative elapsed
+#> 3       analytical    1.000   0.605
+#> 1 deSolve_compiled    1.448   0.876
+#> 2            eigen    2.479   1.500
+# }
+
+# Use stepwise fitting, using optimised parameters from parent only fit, FOMC-SFO
+# \dontrun{
+FOMC_SFO <- mkinmod(
+  parent = mkinsub("FOMC", "m1"),
+  m1 = mkinsub("SFO"))
+#> Temporary DLL for differentials generated and loaded
+fit.FOMC_SFO <- mkinfit(FOMC_SFO, FOCUS_D, quiet = TRUE)
+# Again, we get a warning and try a more sophisticated error model
+fit.FOMC_SFO.tc <- mkinfit(FOMC_SFO, FOCUS_D, quiet = TRUE, error_model = "tc")
+# This model has a higher likelihood, but not significantly so
+lrtest(fit.tc, fit.FOMC_SFO.tc)
+#> Likelihood ratio test
+#> 
+#> Model 1: FOMC_SFO with error model tc and fixed parameter(s) m1_0
+#> Model 2: SFO_SFO with error model tc and fixed parameter(s) m1_0
+#>   #Df  LogLik Df  Chisq Pr(>Chisq)
+#> 1   7 -64.829                     
+#> 2   6 -64.983 -1 0.3075     0.5792
+# Also, the missing standard error for log_beta and the t-tests for alpha
+# and beta indicate overparameterisation
+summary(fit.FOMC_SFO.tc, data = FALSE)
+#> Warning: NaNs produced
+#> Warning: NaNs produced
+#> Warning: diag(.) had 0 or NA entries; non-finite result is doubtful
+#> mkin version used for fitting:    1.1.2 
+#> R version used for fitting:       4.2.1 
+#> Date of fit:     Wed Aug 10 15:24:27 2022 
+#> Date of summary: Wed Aug 10 15:24:27 2022 
+#> 
+#> Equations:
+#> d_parent/dt = - (alpha/beta) * 1/((time/beta) + 1) * parent
+#> d_m1/dt = + f_parent_to_m1 * (alpha/beta) * 1/((time/beta) + 1) *
+#>            parent - k_m1 * m1
+#> 
+#> Model predictions using solution type deSolve 
+#> 
+#> Fitted using 3729 model solutions performed in 2.475 s
+#> 
+#> Error model: Two-component variance function 
+#> 
+#> Error model algorithm: d_3 
+#> Direct fitting and three-step fitting yield approximately the same likelihood 
+#> 
+#> Starting values for parameters to be optimised:
+#>                 value   type
+#> parent_0       100.75  state
+#> alpha            1.00 deparm
+#> beta            10.00 deparm
+#> k_m1             0.10 deparm
+#> f_parent_to_m1   0.50 deparm
+#> sigma_low        0.10  error
+#> rsd_high         0.10  error
+#> 
+#> Starting values for the transformed parameters actually optimised:
+#>                      value lower upper
+#> parent_0        100.750000  -Inf   Inf
+#> log_k_m1         -2.302585  -Inf   Inf
+#> f_parent_qlogis   0.000000  -Inf   Inf
+#> log_alpha         0.000000  -Inf   Inf
+#> log_beta          2.302585  -Inf   Inf
+#> sigma_low         0.100000     0   Inf
+#> rsd_high          0.100000     0   Inf
+#> 
+#> Fixed parameter values:
+#>      value  type
+#> m1_0     0 state
+#> 
+#> Results:
+#> 
+#>       AIC      BIC    logLik
+#>   143.658 155.1211 -64.82902
+#> 
+#> Optimised, transformed parameters with symmetric confidence intervals:
+#>                   Estimate Std. Error     Lower      Upper
+#> parent_0        101.600000  2.6400000 96.240000 107.000000
+#> log_k_m1         -5.284000  0.0929100 -5.474000  -5.095000
+#> f_parent_qlogis   0.001426  0.0767000 -0.155000   0.157800
+#> log_alpha         5.522000  0.0077320  5.506000   5.538000
+#> log_beta          7.806000        NaN       NaN        NaN
+#> sigma_low         0.002488  0.0002431  0.001992   0.002984
+#> rsd_high          0.079210  0.0093280  0.060180   0.098230
+#> 
+#> Parameter correlation:
+#>                  parent_0  log_k_m1 f_parent_qlogis log_alpha log_beta
+#> parent_0         1.000000 -0.095226        -0.76678   0.70544      NaN
+#> log_k_m1        -0.095226  1.000000         0.51432  -0.14387      NaN
+#> f_parent_qlogis -0.766780  0.514321         1.00000  -0.61396      NaN
+#> log_alpha        0.705444 -0.143872        -0.61396   1.00000      NaN
+#> log_beta              NaN       NaN             NaN       NaN        1
+#> sigma_low        0.016073  0.001586         0.01548   5.87007      NaN
+#> rsd_high         0.006626 -0.011700        -0.05357   0.04849      NaN
+#>                 sigma_low  rsd_high
+#> parent_0         0.016073  0.006626
+#> log_k_m1         0.001586 -0.011700
+#> f_parent_qlogis  0.015476 -0.053566
+#> log_alpha        5.870075  0.048487
+#> log_beta              NaN       NaN
+#> sigma_low        1.000000 -0.652558
+#> rsd_high        -0.652558  1.000000
+#> 
+#> Backtransformed parameters:
+#> Confidence intervals for internally transformed parameters are asymmetric.
+#> t-test (unrealistically) based on the assumption of normal distribution
+#> for estimators of untransformed parameters.
+#>                 Estimate t value    Pr(>t)     Lower     Upper
+#> parent_0       1.016e+02 32.7800 6.311e-26 9.624e+01 1.070e+02
+#> k_m1           5.072e-03 10.1200 1.216e-11 4.196e-03 6.130e-03
+#> f_parent_to_m1 5.004e-01 20.8300 4.317e-20 4.613e-01 5.394e-01
+#> alpha          2.502e+02  0.5624 2.889e-01 2.463e+02 2.542e+02
+#> beta           2.455e+03  0.5549 2.915e-01        NA        NA
+#> sigma_low      2.488e-03  0.4843 3.158e-01 1.992e-03 2.984e-03
+#> rsd_high       7.921e-02  8.4300 8.001e-10 6.018e-02 9.823e-02
+#> 
+#> FOCUS Chi2 error levels in percent:
+#>          err.min n.optim df
+#> All data   6.781       5 14
+#> parent     7.141       3  6
+#> m1         4.640       2  8
+#> 
+#> Resulting formation fractions:
+#>                 ff
+#> parent_m1   0.5004
+#> parent_sink 0.4996
+#> 
+#> Estimated disappearance times:
+#>           DT50  DT90 DT50back
+#> parent   6.812  22.7    6.834
+#> m1     136.661 454.0       NA
+
+# We can easily use starting parameters from the parent only fit (only for illustration)
+fit.FOMC = mkinfit("FOMC", FOCUS_2006_D, quiet = TRUE, error_model = "tc")
+fit.FOMC_SFO <- mkinfit(FOMC_SFO, FOCUS_D, quiet = TRUE,
+  parms.ini = fit.FOMC$bparms.ode, error_model = "tc")
+# }
+
    +
    +
    -
    -
    -

    Developed by Johannes Ranke.

    +
    +

    Developed by Johannes Ranke.

    -

    Site built with pkgdown 1.6.1.

    +

    Site built with pkgdown 2.0.6.

    -
    -
    +
    - - + + diff --git a/docs/dev/reference/nlme.mmkin-1.png b/docs/dev/reference/nlme.mmkin-1.png index d93b0950..45e4eebe 100644 Binary files a/docs/dev/reference/nlme.mmkin-1.png and b/docs/dev/reference/nlme.mmkin-1.png differ diff --git a/docs/dev/reference/nlme.mmkin-2.png b/docs/dev/reference/nlme.mmkin-2.png index 5afccfbb..b9a68e92 100644 Binary files a/docs/dev/reference/nlme.mmkin-2.png and b/docs/dev/reference/nlme.mmkin-2.png differ diff --git a/docs/dev/reference/nlme.mmkin-3.png b/docs/dev/reference/nlme.mmkin-3.png index bf3aac6c..2a724bed 100644 Binary files a/docs/dev/reference/nlme.mmkin-3.png and b/docs/dev/reference/nlme.mmkin-3.png differ diff --git a/docs/dev/reference/nlme.mmkin.html b/docs/dev/reference/nlme.mmkin.html index e138ddd4..2bbadb88 100644 --- a/docs/dev/reference/nlme.mmkin.html +++ b/docs/dev/reference/nlme.mmkin.html @@ -19,7 +19,7 @@ have been obtained by fitting the same model to a list of datasets."> mkin - 1.1.0 + 1.1.2 @@ -28,7 +28,7 @@ have been obtained by fitting the same model to a list of datasets.">Functions and data
  • - Example evaluation of FOCUS Laboratory Data L1 to L3
    • +
  • + Example evaluations of dimethenamid data from 2018 with nonlinear mixed-effects models +
  • Example evaluation of FOCUS Example Dataset Z
    • @@ -88,78 +91,115 @@ have been obtained by fitting the same model to a list of datasets.

    - 
    # S3 method for mmkin
-nlme(
-  model,
-  data = "auto",
-  fixed = lapply(as.list(names(mean_degparms(model))), function(el) eval(parse(text =
-    paste(el, 1, sep = "~")))),
-  random = pdDiag(fixed),
-  groups,
-  start = mean_degparms(model, random = TRUE, test_log_parms = TRUE),
-  correlation = NULL,
-  weights = NULL,
-  subset,
-  method = c("ML", "REML"),
-  na.action = na.fail,
-  naPattern,
-  control = list(),
-  verbose = FALSE
-)
-
-# S3 method for nlme.mmkin
-print(x, digits = max(3, getOption("digits") - 3), ...)
-
-# S3 method for nlme.mmkin
-update(object, ...)
    + 
    # S3 method for mmkin
+nlme(
+  model,
+  data = "auto",
+  fixed = lapply(as.list(names(mean_degparms(model))), function(el) eval(parse(text =
+    paste(el, 1, sep = "~")))),
+  random = pdDiag(fixed),
+  groups,
+  start = mean_degparms(model, random = TRUE, test_log_parms = TRUE),
+  correlation = NULL,
+  weights = NULL,
+  subset,
+  method = c("ML", "REML"),
+  na.action = na.fail,
+  naPattern,
+  control = list(),
+  verbose = FALSE
+)
+
+# S3 method for nlme.mmkin
+print(x, digits = max(3, getOption("digits") - 3), ...)
+
+# S3 method for nlme.mmkin
+update(object, ...)

    Arguments

    model

    An mmkin row object.

    + +
    data

    Ignored, data are taken from the mmkin model

    + +
    fixed

    Ignored, all degradation parameters fitted in the mmkin model are used as fixed parameters

    + +
    random

    If not specified, no correlations between random effects are set up for the optimised degradation model parameters. This is achieved by using the nlme::pdDiag method.

    + +
    groups

    See the documentation of nlme

    + +
    start

    If not specified, mean values of the fitted degradation parameters taken from the mmkin object are used

    + +
    correlation

    See the documentation of nlme

    + +
    weights

    passed to nlme

    + +
    subset

    passed to nlme

    + +
    method

    passed to nlme

    + +
    na.action

    passed to nlme

    + +
    naPattern

    passed to nlme

    + +
    control

    passed to nlme

    + +
    verbose

    passed to nlme

    + +
    x

    An nlme.mmkin object to print

    + +
    digits

    Number of digits to use for printing

    + +
    ...

    Update specifications passed to update.nlme

    + +
    object

    An nlme.mmkin object to update

    +

    Value

    -

    Upon success, a fitted 'nlme.mmkin' object, which is an nlme object + + +

    Upon success, a fitted 'nlme.mmkin' object, which is an nlme object with additional elements. It also inherits from 'mixed.mmkin'.

    @@ -182,20 +222,19 @@ methods that will automatically work on 'nlme.mmkin' objects, such as

    Examples

    - 
    ds <- lapply(experimental_data_for_UBA_2019[6:10],
- function(x) subset(x$data[c("name", "time", "value")], name == "parent"))
-f <- mmkin(c("SFO", "DFOP"), ds, quiet = TRUE, cores = 1)
-library(nlme)
-f_nlme_sfo <- nlme(f["SFO", ])
-
-# \dontrun{
-
-  f_nlme_dfop <- nlme(f["DFOP", ])
-  anova(f_nlme_sfo, f_nlme_dfop)
+    
    ds <- lapply(experimental_data_for_UBA_2019[6:10],
+ function(x) subset(x$data[c("name", "time", "value")], name == "parent"))
+
+# \dontrun{
+  f <- mmkin(c("SFO", "DFOP"), ds, quiet = TRUE, cores = 1)
+  library(nlme)
+  f_nlme_sfo <- nlme(f["SFO", ])
+  f_nlme_dfop <- nlme(f["DFOP", ])
+  anova(f_nlme_sfo, f_nlme_dfop)
 #>             Model df      AIC      BIC    logLik   Test  L.Ratio p-value
 #> f_nlme_sfo      1  5 625.0539 637.5529 -307.5269                        
 #> f_nlme_dfop     2  9 495.1270 517.6253 -238.5635 1 vs 2 137.9269  <.0001
-  print(f_nlme_dfop)
+  print(f_nlme_dfop)
 #> Kinetic nonlinear mixed-effects model fit by maximum likelihood
 #> 
 #> Structural model:
@@ -220,50 +259,50 @@ methods that will automatically work on 'nlme.mmkin' objects, such as
 #>         parent_0 log_k1 log_k2 g_qlogis Residual
 #> StdDev:    2.488 0.8447   1.33   0.4652    2.321
 #> 
-  plot(f_nlme_dfop)
+  plot(f_nlme_dfop)
 
-  endpoints(f_nlme_dfop)
+  endpoints(f_nlme_dfop)
 #> $distimes
 #>            DT50     DT90 DT50back  DT50_k1  DT50_k2
 #> parent 10.79857 100.7937 30.34193 4.193938 43.85443
 #> 
-
-  ds_2 <- lapply(experimental_data_for_UBA_2019[6:10],
-   function(x) x$data[c("name", "time", "value")])
-  m_sfo_sfo <- mkinmod(parent = mkinsub("SFO", "A1"),
-    A1 = mkinsub("SFO"), use_of_ff = "min", quiet = TRUE)
-  m_sfo_sfo_ff <- mkinmod(parent = mkinsub("SFO", "A1"),
-    A1 = mkinsub("SFO"), use_of_ff = "max", quiet = TRUE)
-  m_dfop_sfo <- mkinmod(parent = mkinsub("DFOP", "A1"),
-    A1 = mkinsub("SFO"), quiet = TRUE)
-
-  f_2 <- mmkin(list("SFO-SFO" = m_sfo_sfo,
-   "SFO-SFO-ff" = m_sfo_sfo_ff,
-   "DFOP-SFO" = m_dfop_sfo),
-    ds_2, quiet = TRUE)
-
-  f_nlme_sfo_sfo <- nlme(f_2["SFO-SFO", ])
-  plot(f_nlme_sfo_sfo)
+
+  ds_2 <- lapply(experimental_data_for_UBA_2019[6:10],
+   function(x) x$data[c("name", "time", "value")])
+  m_sfo_sfo <- mkinmod(parent = mkinsub("SFO", "A1"),
+    A1 = mkinsub("SFO"), use_of_ff = "min", quiet = TRUE)
+  m_sfo_sfo_ff <- mkinmod(parent = mkinsub("SFO", "A1"),
+    A1 = mkinsub("SFO"), use_of_ff = "max", quiet = TRUE)
+  m_dfop_sfo <- mkinmod(parent = mkinsub("DFOP", "A1"),
+    A1 = mkinsub("SFO"), quiet = TRUE)
+
+  f_2 <- mmkin(list("SFO-SFO" = m_sfo_sfo,
+   "SFO-SFO-ff" = m_sfo_sfo_ff,
+   "DFOP-SFO" = m_dfop_sfo),
+    ds_2, quiet = TRUE)
+
+  f_nlme_sfo_sfo <- nlme(f_2["SFO-SFO", ])
+  plot(f_nlme_sfo_sfo)
 
-
-  # With formation fractions this does not coverge with defaults
-  # f_nlme_sfo_sfo_ff <- nlme(f_2["SFO-SFO-ff", ])
-  #plot(f_nlme_sfo_sfo_ff)
-
-  # For the following, we need to increase pnlsMaxIter and the tolerance
-  # to get convergence
-  f_nlme_dfop_sfo <- nlme(f_2["DFOP-SFO", ],
-    control = list(pnlsMaxIter = 120, tolerance = 5e-4))
-
-  plot(f_nlme_dfop_sfo)
+
+  # With formation fractions this does not coverge with defaults
+  # f_nlme_sfo_sfo_ff <- nlme(f_2["SFO-SFO-ff", ])
+  #plot(f_nlme_sfo_sfo_ff)
+
+  # For the following, we need to increase pnlsMaxIter and the tolerance
+  # to get convergence
+  f_nlme_dfop_sfo <- nlme(f_2["DFOP-SFO", ],
+    control = list(pnlsMaxIter = 120, tolerance = 5e-4))
+
+  plot(f_nlme_dfop_sfo)
 
-
-  anova(f_nlme_dfop_sfo, f_nlme_sfo_sfo)
+
+  anova(f_nlme_dfop_sfo, f_nlme_sfo_sfo)
 #>                 Model df       AIC       BIC    logLik   Test  L.Ratio p-value
 #> f_nlme_dfop_sfo     1 13  843.8547  884.6201 -408.9274                        
 #> f_nlme_sfo_sfo      2  9 1085.1821 1113.4043 -533.5910 1 vs 2 249.3274  <.0001
-
-  endpoints(f_nlme_sfo_sfo)
+
+  endpoints(f_nlme_sfo_sfo)
 #> $ff
 #> parent_sink   parent_A1     A1_sink 
 #>   0.5912432   0.4087568   1.0000000 
@@ -273,7 +312,7 @@ methods that will automatically work on 'nlme.mmkin' objects, such as
 #> parent 19.13518  63.5657
 #> A1     66.02155 219.3189
 #> 
-  endpoints(f_nlme_dfop_sfo)
+  endpoints(f_nlme_dfop_sfo)
 #> $ff
 #>   parent_A1 parent_sink 
 #>   0.2768574   0.7231426 
@@ -283,17 +322,17 @@ methods that will automatically work on 'nlme.mmkin' objects, such as
 #> parent  11.07091 104.6320 31.49737 4.462383 46.20825
 #> A1     162.30519 539.1662       NA       NA       NA
 #> 
-
-  if (length(findFunction("varConstProp")) > 0) { # tc error model for nlme available
-    # Attempts to fit metabolite kinetics with the tc error model are possible,
-    # but need tweeking of control values and sometimes do not converge
-
-    f_tc <- mmkin(c("SFO", "DFOP"), ds, quiet = TRUE, error_model = "tc")
-    f_nlme_sfo_tc <- nlme(f_tc["SFO", ])
-    f_nlme_dfop_tc <- nlme(f_tc["DFOP", ])
-    AIC(f_nlme_sfo, f_nlme_sfo_tc, f_nlme_dfop, f_nlme_dfop_tc)
-    print(f_nlme_dfop_tc)
-  }
+
+  if (length(findFunction("varConstProp")) > 0) { # tc error model for nlme available
+    # Attempts to fit metabolite kinetics with the tc error model are possible,
+    # but need tweeking of control values and sometimes do not converge
+
+    f_tc <- mmkin(c("SFO", "DFOP"), ds, quiet = TRUE, error_model = "tc")
+    f_nlme_sfo_tc <- nlme(f_tc["SFO", ])
+    f_nlme_dfop_tc <- nlme(f_tc["DFOP", ])
+    AIC(f_nlme_sfo, f_nlme_sfo_tc, f_nlme_dfop, f_nlme_dfop_tc)
+    print(f_nlme_dfop_tc)
+  }
 #> Kinetic nonlinear mixed-effects model fit by maximum likelihood
 #> 
 #> Structural model:
@@ -324,10 +363,10 @@ methods that will automatically work on 'nlme.mmkin' objects, such as
 #>  Parameter estimates:
 #>      const       prop 
 #> 2.23223147 0.01262395 
-
-  f_2_obs <- update(f_2, error_model = "obs")
-  f_nlme_sfo_sfo_obs <- nlme(f_2_obs["SFO-SFO", ])
-  print(f_nlme_sfo_sfo_obs)
+
+  f_2_obs <- update(f_2, error_model = "obs")
+  f_nlme_sfo_sfo_obs <- nlme(f_2_obs["SFO-SFO", ])
+  print(f_nlme_sfo_sfo_obs)
 #> Kinetic nonlinear mixed-effects model fit by maximum likelihood
 #> 
 #> Structural model:
@@ -357,23 +396,23 @@ methods that will automatically work on 'nlme.mmkin' objects, such as
 #>  Parameter estimates:
 #>    parent        A1 
 #> 1.0000000 0.2049995 
-  f_nlme_dfop_sfo_obs <- nlme(f_2_obs["DFOP-SFO", ],
-    control = list(pnlsMaxIter = 120, tolerance = 5e-4))
-
-  f_2_tc <- update(f_2, error_model = "tc")
-  # f_nlme_sfo_sfo_tc <- nlme(f_2_tc["SFO-SFO", ]) # No convergence with 50 iterations
-  # f_nlme_dfop_sfo_tc <- nlme(f_2_tc["DFOP-SFO", ],
-  #  control = list(pnlsMaxIter = 120, tolerance = 5e-4)) # Error in X[, fmap[[nm]]] <- gradnm
-
-  anova(f_nlme_dfop_sfo, f_nlme_dfop_sfo_obs)
+  f_nlme_dfop_sfo_obs <- nlme(f_2_obs["DFOP-SFO", ],
+    control = list(pnlsMaxIter = 120, tolerance = 5e-4))
+
+  f_2_tc <- update(f_2, error_model = "tc")
+  # f_nlme_sfo_sfo_tc <- nlme(f_2_tc["SFO-SFO", ]) # No convergence with 50 iterations
+  # f_nlme_dfop_sfo_tc <- nlme(f_2_tc["DFOP-SFO", ],
+  #  control = list(pnlsMaxIter = 120, tolerance = 5e-4)) # Error in X[, fmap[[nm]]] <- gradnm
+
+  anova(f_nlme_dfop_sfo, f_nlme_dfop_sfo_obs)
 #>                     Model df      AIC      BIC    logLik   Test  L.Ratio
 #> f_nlme_dfop_sfo         1 13 843.8547 884.6201 -408.9274                
 #> f_nlme_dfop_sfo_obs     2 14 817.5338 861.4350 -394.7669 1 vs 2 28.32091
 #>                     p-value
 #> f_nlme_dfop_sfo            
 #> f_nlme_dfop_sfo_obs  <.0001
-
-# }
+
+# }
    @@ -388,7 +427,7 @@ methods that will automatically work on 'nlme.mmkin' objects, such as
    -

    Site built with pkgdown 2.0.2.

    +

    Site built with pkgdown 2.0.6.

    diff --git a/docs/dev/reference/plot.mixed.mmkin-2.png b/docs/dev/reference/plot.mixed.mmkin-2.png index 67c9b8ab..8678c166 100644 Binary files a/docs/dev/reference/plot.mixed.mmkin-2.png and b/docs/dev/reference/plot.mixed.mmkin-2.png differ diff --git a/docs/dev/reference/plot.mixed.mmkin-3.png b/docs/dev/reference/plot.mixed.mmkin-3.png index 7e2876b3..9bd01852 100644 Binary files a/docs/dev/reference/plot.mixed.mmkin-3.png and b/docs/dev/reference/plot.mixed.mmkin-3.png differ diff --git a/docs/dev/reference/plot.mixed.mmkin-4.png b/docs/dev/reference/plot.mixed.mmkin-4.png index 945c4d41..d2fa1ae1 100644 Binary files a/docs/dev/reference/plot.mixed.mmkin-4.png and b/docs/dev/reference/plot.mixed.mmkin-4.png differ diff --git a/docs/dev/reference/plot.mixed.mmkin.html b/docs/dev/reference/plot.mixed.mmkin.html index 19b8e2c1..83017908 100644 --- a/docs/dev/reference/plot.mixed.mmkin.html +++ b/docs/dev/reference/plot.mixed.mmkin.html @@ -17,7 +17,7 @@ mkin - 1.1.0 + 1.1.2 @@ -26,7 +26,7 @@ Functions and data
  • -
  • Example evaluation of FOCUS Laboratory Data L1 to L3
    • +
  • + Example evaluations of dimethenamid data from 2018 with nonlinear mixed-effects models +
  • Example evaluation of FOCUS Example Dataset Z
    • @@ -84,90 +87,139 @@
    - 
    # S3 method for mixed.mmkin
-plot(
-  x,
-  i = 1:ncol(x$mmkin),
-  obs_vars = names(x$mkinmod$map),
-  standardized = TRUE,
-  xlab = "Time",
-  xlim = range(x$data$time),
-  resplot = c("predicted", "time"),
-  pred_over = NULL,
-  test_log_parms = FALSE,
-  conf.level = 0.6,
-  ymax = "auto",
-  maxabs = "auto",
-  ncol.legend = ifelse(length(i) <= 3, length(i) + 1, ifelse(length(i) <= 8, 3, 4)),
-  nrow.legend = ceiling((length(i) + 1)/ncol.legend),
-  rel.height.legend = 0.02 + 0.07 * nrow.legend,
-  rel.height.bottom = 1.1,
-  pch_ds = 1:length(i),
-  col_ds = pch_ds + 1,
-  lty_ds = col_ds,
-  frame = TRUE,
-  ...
-)
    + 
    # S3 method for mixed.mmkin
+plot(
+  x,
+  i = 1:ncol(x$mmkin),
+  obs_vars = names(x$mkinmod$map),
+  standardized = TRUE,
+  xlab = "Time",
+  xlim = range(x$data$time),
+  resplot = c("predicted", "time"),
+  pred_over = NULL,
+  test_log_parms = FALSE,
+  conf.level = 0.6,
+  default_log_parms = NA,
+  ymax = "auto",
+  maxabs = "auto",
+  ncol.legend = ifelse(length(i) <= 3, length(i) + 1, ifelse(length(i) <= 8, 3, 4)),
+  nrow.legend = ceiling((length(i) + 1)/ncol.legend),
+  rel.height.legend = 0.02 + 0.07 * nrow.legend,
+  rel.height.bottom = 1.1,
+  pch_ds = 1:length(i),
+  col_ds = pch_ds + 1,
+  lty_ds = col_ds,
+  frame = TRUE,
+  ...
+)

    Arguments

    x

    An object of class mixed.mmkin, saem.mmkin or nlme.mmkin

    + +
    i

    A numeric index to select datasets for which to plot the individual predictions, in case plots get too large

    + +
    obs_vars

    A character vector of names of the observed variables for which the data and the model should be plotted. Defauls to all observed variables in the model.

    + +
    standardized

    Should the residuals be standardized? Only takes effect if resplot = "time".

    + +
    xlab

    Label for the x axis.

    + +
    xlim

    Plot range in x direction.

    + +
    resplot

    Should the residuals plotted against time or against predicted values?

    + +
    pred_over

    Named list of alternative predictions as obtained from mkinpredict with a compatible mkinmod.

    + +
    test_log_parms

    Passed to mean_degparms in the case of an mixed.mmkin object

    + +
    conf.level

    Passed to mean_degparms in the case of an mixed.mmkin object

    + + +
    default_log_parms
    +

    Passed to mean_degparms in the case of an +mixed.mmkin object

    + +
    ymax

    Vector of maximum y axis values

    + +
    maxabs

    Maximum absolute value of the residuals. This is used for the scaling of the y axis and defaults to "auto".

    + +
    ncol.legend

    Number of columns to use in the legend

    + +
    nrow.legend

    Number of rows to use in the legend

    + +
    rel.height.legend

    The relative height of the legend shown on top

    + +
    rel.height.bottom

    The relative height of the bottom plot row

    + +
    pch_ds

    Symbols to be used for plotting the data.

    + +
    col_ds

    Colors used for plotting the observed data and the corresponding model prediction lines for the different datasets.

    + +
    lty_ds

    Line types to be used for the model predictions.

    + +
    frame

    Should a frame be drawn around the plots?

    + +
    ...

    Further arguments passed to plot.

    +

    Value

    -

    The function is called for its side effect.

    + + +

    The function is called for its side effect.

    Author

    @@ -176,42 +228,41 @@ corresponding model prediction lines for the different datasets.

    Examples

    - 
    ds <- lapply(experimental_data_for_UBA_2019[6:10],
- function(x) x$data[c("name", "time", "value")])
-names(ds) <- paste0("ds ", 6:10)
-dfop_sfo <- mkinmod(parent = mkinsub("DFOP", "A1"),
-  A1 = mkinsub("SFO"), quiet = TRUE)
-# \dontrun{
-f <- mmkin(list("DFOP-SFO" = dfop_sfo), ds, quiet = TRUE)
-plot(f[, 3:4], standardized = TRUE)
+    
    ds <- lapply(experimental_data_for_UBA_2019[6:10],
+ function(x) x$data[c("name", "time", "value")])
+names(ds) <- paste0("ds ", 6:10)
+dfop_sfo <- mkinmod(parent = mkinsub("DFOP", "A1"),
+  A1 = mkinsub("SFO"), quiet = TRUE)
+# \dontrun{
+f <- mmkin(list("DFOP-SFO" = dfop_sfo), ds, quiet = TRUE)
+plot(f[, 3:4], standardized = TRUE)
 
-
-# For this fit we need to increase pnlsMaxiter, and we increase the
-# tolerance in order to speed up the fit for this example evaluation
-# It still takes 20 seconds to run
-f_nlme <- nlme(f, control = list(pnlsMaxIter = 120, tolerance = 1e-3))
-plot(f_nlme)
+
+# For this fit we need to increase pnlsMaxiter, and we increase the
+# tolerance in order to speed up the fit for this example evaluation
+# It still takes 20 seconds to run
+f_nlme <- nlme(f, control = list(pnlsMaxIter = 120, tolerance = 1e-3))
+plot(f_nlme)
 
-
-f_saem <- saem(f, transformations = "saemix")
-#> Error in saem(f, transformations = "saemix"): unused argument (transformations = "saemix")
-plot(f_saem)
-#> Error in plot(f_saem): object 'f_saem' not found
-
-f_obs <- mmkin(list("DFOP-SFO" = dfop_sfo), ds, quiet = TRUE, error_model = "obs")
-f_nlmix <- nlmix(f_obs)
+
+f_saem <- saem(f, transformations = "saemix")
+plot(f_saem)
+
+
+f_obs <- mmkin(list("DFOP-SFO" = dfop_sfo), ds, quiet = TRUE, error_model = "obs")
+f_nlmix <- nlmix(f_obs)
 #> Error in nlmix(f_obs): could not find function "nlmix"
-plot(f_nlmix)
+plot(f_nlmix)
 #> Error in plot(f_nlmix): object 'f_nlmix' not found
-
-# 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))
-#> Error in plot(f_saem, pred_over = list(nlme = pred_nlme)): object 'f_saem' not found
-# }
+
+# 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))
+
+# }
    @@ -226,7 +277,7 @@ corresponding model prediction lines for the different datasets.

    -

    Site built with pkgdown 2.0.2.

    +

    Site built with pkgdown 2.0.6.

    diff --git a/docs/dev/reference/plot.mkinfit-1.png b/docs/dev/reference/plot.mkinfit-1.png index e5da9f1c..ea1032fb 100644 Binary files a/docs/dev/reference/plot.mkinfit-1.png and b/docs/dev/reference/plot.mkinfit-1.png differ diff --git a/docs/dev/reference/plot.mkinfit-2.png b/docs/dev/reference/plot.mkinfit-2.png index a11d1680..cef94cb8 100644 Binary files a/docs/dev/reference/plot.mkinfit-2.png and b/docs/dev/reference/plot.mkinfit-2.png differ diff --git a/docs/dev/reference/plot.mkinfit-3.png b/docs/dev/reference/plot.mkinfit-3.png index c976d4b1..8a9dbd13 100644 Binary files a/docs/dev/reference/plot.mkinfit-3.png and b/docs/dev/reference/plot.mkinfit-3.png differ diff --git a/docs/dev/reference/plot.mkinfit-4.png b/docs/dev/reference/plot.mkinfit-4.png index c8bc00fe..a7164caa 100644 Binary files a/docs/dev/reference/plot.mkinfit-4.png and b/docs/dev/reference/plot.mkinfit-4.png differ diff --git a/docs/dev/reference/plot.mkinfit-5.png b/docs/dev/reference/plot.mkinfit-5.png index 6631aa68..f90b3f54 100644 Binary files a/docs/dev/reference/plot.mkinfit-5.png and b/docs/dev/reference/plot.mkinfit-5.png differ diff --git a/docs/dev/reference/plot.mkinfit-6.png b/docs/dev/reference/plot.mkinfit-6.png index 946b20c5..3d0fb25e 100644 Binary files a/docs/dev/reference/plot.mkinfit-6.png and b/docs/dev/reference/plot.mkinfit-6.png differ diff --git a/docs/dev/reference/plot.mkinfit-7.png b/docs/dev/reference/plot.mkinfit-7.png index 10807ea8..3e5d828e 100644 Binary files a/docs/dev/reference/plot.mkinfit-7.png and b/docs/dev/reference/plot.mkinfit-7.png differ diff --git a/docs/dev/reference/plot.mkinfit.html b/docs/dev/reference/plot.mkinfit.html index c5249ecd..764f0699 100644 --- a/docs/dev/reference/plot.mkinfit.html +++ b/docs/dev/reference/plot.mkinfit.html @@ -1,69 +1,14 @@ - - - - - - - -Plot the observed data and the fitted model of an mkinfit object — plot.mkinfit • mkin - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -Plot the observed data and the fitted model of an mkinfit object — plot.mkinfit • mkin - - - - - - - - - - - + + - - -
    -
    - -
    - -
    +

    Plot the observed data and the fitted model of an mkinfit object

    - Source: R/plot.mkinfit.R + Source: R/plot.mkinfit.R

    Solves the differential equations with the optimised and fixed parameters -from a previous successful call to mkinfit and plots the +from a previous successful call to mkinfit and plots the observed data together with the solution of the fitted model.

    - 
    # S3 method for mkinfit
-plot(
-  x,
-  fit = x,
-  obs_vars = names(fit$mkinmod$map),
-  xlab = "Time",
-  ylab = "Residue",
-  xlim = range(fit$data$time),
-  ylim = "default",
-  col_obs = 1:length(obs_vars),
-  pch_obs = col_obs,
-  lty_obs = rep(1, length(obs_vars)),
-  add = FALSE,
-  legend = !add,
-  show_residuals = FALSE,
-  show_errplot = FALSE,
-  maxabs = "auto",
-  sep_obs = FALSE,
-  rel.height.middle = 0.9,
-  row_layout = FALSE,
-  lpos = "topright",
-  inset = c(0.05, 0.05),
-  show_errmin = FALSE,
-  errmin_digits = 3,
-  frame = TRUE,
-  ...
-)
-
-plot_sep(
-  fit,
-  show_errmin = TRUE,
-  show_residuals = ifelse(identical(fit$err_mod, "const"), TRUE, "standardized"),
-  ...
-)
-
-plot_res(
-  fit,
-  sep_obs = FALSE,
-  show_errmin = sep_obs,
-  standardized = ifelse(identical(fit$err_mod, "const"), FALSE, TRUE),
-  ...
-)
-
-plot_err(fit, sep_obs = FALSE, show_errmin = sep_obs, ...)
    - -

    Arguments

    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    x

    Alias for fit introduced for compatibility with the generic S3 -method.

    fit

    An object of class mkinfit.

    obs_vars

    A character vector of names of the observed variables for +

    + 
    # S3 method for mkinfit
+plot(
+  x,
+  fit = x,
+  obs_vars = names(fit$mkinmod$map),
+  xlab = "Time",
+  ylab = "Residue",
+  xlim = range(fit$data$time),
+  ylim = "default",
+  col_obs = 1:length(obs_vars),
+  pch_obs = col_obs,
+  lty_obs = rep(1, length(obs_vars)),
+  add = FALSE,
+  legend = !add,
+  show_residuals = FALSE,
+  show_errplot = FALSE,
+  maxabs = "auto",
+  sep_obs = FALSE,
+  rel.height.middle = 0.9,
+  row_layout = FALSE,
+  lpos = "topright",
+  inset = c(0.05, 0.05),
+  show_errmin = FALSE,
+  errmin_digits = 3,
+  frame = TRUE,
+  ...
+)
+
+plot_sep(
+  fit,
+  show_errmin = TRUE,
+  show_residuals = ifelse(identical(fit$err_mod, "const"), TRUE, "standardized"),
+  ...
+)
+
+plot_res(
+  fit,
+  sep_obs = FALSE,
+  show_errmin = sep_obs,
+  standardized = ifelse(identical(fit$err_mod, "const"), FALSE, TRUE),
+  ...
+)
+
+plot_err(fit, sep_obs = FALSE, show_errmin = sep_obs, ...)
    +
    + +
    +

    Arguments

    +
    x
    +

    Alias for fit introduced for compatibility with the generic S3 +method.

    + + +
    fit
    +

    An object of class mkinfit.

    + + +
    obs_vars
    +

    A character vector of names of the observed variables for which the data and the model should be plotted. Defauls to all observed -variables in the model.

    xlab

    Label for the x axis.

    ylab

    Label for the y axis.

    xlim

    Plot range in x direction.

    ylim

    Plot range in y direction.

    col_obs

    Colors used for plotting the observed data and the -corresponding model prediction lines.

    pch_obs

    Symbols to be used for plotting the data.

    lty_obs

    Line types to be used for the model predictions.

    add

    Should the plot be added to an existing plot?

    legend

    Should a legend be added to the plot?

    show_residuals

    Should residuals be shown? If only one plot of the +variables in the model.

    + + +
    xlab
    +

    Label for the x axis.

    + + +
    ylab
    +

    Label for the y axis.

    + + +
    xlim
    +

    Plot range in x direction.

    + + +
    ylim
    +

    Plot range in y direction. If given as a list, plot ranges +for the different plot rows can be given for row layout.

    + + +
    col_obs
    +

    Colors used for plotting the observed data and the +corresponding model prediction lines.

    + + +
    pch_obs
    +

    Symbols to be used for plotting the data.

    + + +
    lty_obs
    +

    Line types to be used for the model predictions.

    + + +
    add
    +

    Should the plot be added to an existing plot?

    + + +
    legend
    +

    Should a legend be added to the plot?

    + + +
    show_residuals
    +

    Should residuals be shown? If only one plot of the fits is shown, the residual plot is in the lower third of the plot. Otherwise, i.e. if "sep_obs" is given, the residual plots will be located to the right of the plots of the fitted curves. If this is set to 'standardized', a plot of the residuals divided by the standard deviation -given by the fitted error model will be shown.

    show_errplot

    Should squared residuals and the error model be shown? +given by the fitted error model will be shown.

    + + +
    show_errplot
    +

    Should squared residuals and the error model be shown? If only one plot of the fits is shown, this plot is in the lower third of the plot. Otherwise, i.e. if "sep_obs" is given, the residual plots will -be located to the right of the plots of the fitted curves.

    maxabs

    Maximum absolute value of the residuals. This is used for the -scaling of the y axis and defaults to "auto".

    sep_obs

    Should the observed variables be shown in separate subplots? +be located to the right of the plots of the fitted curves.

    + + +
    maxabs
    +

    Maximum absolute value of the residuals. This is used for the +scaling of the y axis and defaults to "auto".

    + + +
    sep_obs
    +

    Should the observed variables be shown in separate subplots? If yes, residual plots requested by "show_residuals" will be shown next -to, not below the plot of the fits.

    rel.height.middle

    The relative height of the middle plot, if more -than two rows of plots are shown.

    row_layout

    Should we use a row layout where the residual plot or the -error model plot is shown to the right?

    lpos

    Position(s) of the legend(s). Passed to legend as +to, not below the plot of the fits.

    + + +
    rel.height.middle
    +

    The relative height of the middle plot, if more +than two rows of plots are shown.

    + + +
    row_layout
    +

    Should we use a row layout where the residual plot or the +error model plot is shown to the right?

    + + +
    lpos
    +

    Position(s) of the legend(s). Passed to legend as the first argument. If not length one, this should be of the same length -as the obs_var argument.

    inset

    Passed to legend if applicable.

    show_errmin

    Should the FOCUS chi2 error value be shown in the upper -margin of the plot?

    errmin_digits

    The number of significant digits for rounding the FOCUS -chi2 error percentage.

    frame

    Should a frame be drawn around the plots?

    ...

    Further arguments passed to plot.

    standardized

    When calling 'plot_res', should the residuals be -standardized in the residual plot?

    - -

    Value

    - -

    The function is called for its side effect.

    -

    Details

    - -

    If the current plot device is a tikz device, then +as the obs_var argument.

    + + +
    inset
    +

    Passed to legend if applicable.

    + + +
    show_errmin
    +

    Should the FOCUS chi2 error value be shown in the upper +margin of the plot?

    + + +
    errmin_digits
    +

    The number of significant digits for rounding the FOCUS +chi2 error percentage.

    + + +
    frame
    +

    Should a frame be drawn around the plots?

    + + +
    ...
    +

    Further arguments passed to plot.

    + + +
    standardized
    +

    When calling 'plot_res', should the residuals be +standardized in the residual plot?

    + +
    +
    +

    Value

    + + +

    The function is called for its side effect.

    +
    +
    +

    Details

    +

    If the current plot device is a tikz device, then latex is being used for the formatting of the chi2 error level, if show_errmin = TRUE.

    -

    Author

    - +
    +
    +

    Author

    Johannes Ranke

    +
    -

    Examples

    - 
    
-# One parent compound, one metabolite, both single first order, path from
-# parent to sink included
-# \dontrun{
-SFO_SFO <- mkinmod(parent = mkinsub("SFO", "m1", full = "Parent"),
-                   m1 = mkinsub("SFO", full = "Metabolite M1" ))
-
    #> Temporary DLL for differentials generated and loaded
    fit <- mkinfit(SFO_SFO, FOCUS_2006_D, quiet = TRUE)
-
    #> Warning: Observations with value of zero were removed from the data
    fit <- mkinfit(SFO_SFO, FOCUS_2006_D, quiet = TRUE, error_model = "tc")
-
    #> Warning: Observations with value of zero were removed from the data
    plot(fit)
-
    plot_res(fit)
-
    plot_res(fit, standardized = FALSE)
-
    plot_err(fit)
-
    
-# Show the observed variables separately, with residuals
-plot(fit, sep_obs = TRUE, show_residuals = TRUE, lpos = c("topright", "bottomright"),
-     show_errmin = TRUE)
-
    
-# The same can be obtained with less typing, using the convenience function plot_sep
-plot_sep(fit, lpos = c("topright", "bottomright"))
-
    
-# Show the observed variables separately, with the error model
-plot(fit, sep_obs = TRUE, show_errplot = TRUE, lpos = c("topright", "bottomright"),
-     show_errmin = TRUE)
-
    # }
-
-
    +
    +

    Examples

    + 
    
+# One parent compound, one metabolite, both single first order, path from
+# parent to sink included
+# \dontrun{
+SFO_SFO <- mkinmod(parent = mkinsub("SFO", "m1", full = "Parent"),
+                   m1 = mkinsub("SFO", full = "Metabolite M1" ))
+#> Temporary DLL for differentials generated and loaded
+fit <- mkinfit(SFO_SFO, FOCUS_2006_D, quiet = TRUE)
+#> Warning: Observations with value of zero were removed from the data
+fit <- mkinfit(SFO_SFO, FOCUS_2006_D, quiet = TRUE, error_model = "tc")
+#> Warning: Observations with value of zero were removed from the data
+plot(fit)
+
+plot_res(fit)
+
+plot_res(fit, standardized = FALSE)
+
+plot_err(fit)
+
+
+# Show the observed variables separately, with residuals
+plot(fit, sep_obs = TRUE, show_residuals = TRUE, lpos = c("topright", "bottomright"),
+     show_errmin = TRUE)
+
+
+# The same can be obtained with less typing, using the convenience function plot_sep
+plot_sep(fit, lpos = c("topright", "bottomright"))
+
+
+# Show the observed variables separately, with the error model
+plot(fit, sep_obs = TRUE, show_errplot = TRUE, lpos = c("topright", "bottomright"),
+     show_errmin = TRUE)
+
+# }
+
+
    +
    +
    -
    -
    -

    Developed by Johannes Ranke.

    +
    +

    Developed by Johannes Ranke.

    -

    Site built with pkgdown 1.6.1.

    +

    Site built with pkgdown 2.0.6.

    -
    -
    +
    - - + + diff --git a/docs/dev/reference/reexports.html b/docs/dev/reference/reexports.html index b8ead4b5..0999e346 100644 --- a/docs/dev/reference/reexports.html +++ b/docs/dev/reference/reexports.html @@ -10,10 +10,6 @@ lrtest intervals, nlme - nlmixr -nlmixr - - "> - - - - - - -Summary method for class "mkinfit" — summary.mkinfit • mkin - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -Summary method for class "mkinfit" — summary.mkinfit • mkin - - - - - - - - - - - + + - - -
    -
    - -
    - -
    +

    Summary method for class "mkinfit"

    - Source: R/summary.mkinfit.R + Source: R/summary.mkinfit.R
    @@ -155,187 +94,217 @@ and optionally the data, consisting of observed, predicted and residual values.

    - 
    # S3 method for mkinfit
-summary(object, data = TRUE, distimes = TRUE, alpha = 0.05, ...)
-
-# S3 method for summary.mkinfit
-print(x, digits = max(3, getOption("digits") - 3), ...)
    - -

    Arguments

    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    object

    an object of class mkinfit.

    data

    logical, indicating whether the data should be included in the -summary.

    distimes

    logical, indicating whether DT50 and DT90 values should be -included.

    alpha

    error level for confidence interval estimation from t -distribution

    ...

    optional arguments passed to methods like print.

    x

    an object of class summary.mkinfit.

    digits

    Number of digits to use for printing

    - -

    Value

    - -

    The summary function returns a list with components, among others

    -
    version, Rversion

    The mkin and R versions used

    -
    date.fit, date.summary

    The dates where the fit and the summary were +

    + 
    # S3 method for mkinfit
+summary(object, data = TRUE, distimes = TRUE, alpha = 0.05, ...)
+
+# S3 method for summary.mkinfit
+print(x, digits = max(3, getOption("digits") - 3), ...)
    +
    + +
    +

    Arguments

    +
    object
    +

    an object of class mkinfit.

    + + +
    data
    +

    logical, indicating whether the data should be included in the +summary.

    + + +
    distimes
    +

    logical, indicating whether DT50 and DT90 values should be +included.

    + + +
    alpha
    +

    error level for confidence interval estimation from t +distribution

    + + +
    ...
    +

    optional arguments passed to methods like print.

    + + +
    x
    +

    an object of class summary.mkinfit.

    + + +
    digits
    +

    Number of digits to use for printing

    + +
    +
    +

    Value

    + + +

    The summary function returns a list with components, among others

    +
    version, Rversion
    +

    The mkin and R versions used

    + +
    date.fit, date.summary
    +

    The dates where the fit and the summary were produced

    -
    diffs

    The differential equations used in the model

    -
    use_of_ff

    Was maximum or minimum use made of formation fractions

    -
    bpar

    Optimised and backtransformed + +

    diffs
    +

    The differential equations used in the model

    + +
    use_of_ff
    +

    Was maximum or minimum use made of formation fractions

    + +
    bpar
    +

    Optimised and backtransformed parameters

    -
    data

    The data (see Description above).

    -
    start

    The starting values and bounds, if applicable, for optimised + +

    data
    +

    The data (see Description above).

    + +
    start
    +

    The starting values and bounds, if applicable, for optimised parameters.

    -
    fixed

    The values of fixed parameters.

    -
    errmin

    The chi2 error levels for + +

    fixed
    +

    The values of fixed parameters.

    + +
    errmin
    +

    The chi2 error levels for each observed variable.

    -
    bparms.ode

    All backtransformed ODE + +

    bparms.ode
    +

    All backtransformed ODE parameters, for use as starting parameters for related models.

    -
    errparms

    Error model parameters.

    -
    ff

    The estimated formation fractions derived from the fitted + +

    errparms
    +

    Error model parameters.

    + +
    ff
    +

    The estimated formation fractions derived from the fitted model.

    -
    distimes

    The DT50 and DT90 values for each observed variable.

    -
    SFORB

    If applicable, eigenvalues of SFORB components of the model.

    -The print method is called for its side effect, i.e. printing the summary. -

    References

    +
    distimes
    +

    The DT50 and DT90 values for each observed variable.

    + +
    SFORB
    +

    If applicable, eigenvalues of SFORB components of the model.

    -

    FOCUS (2006) “Guidance Document on Estimating Persistence +

    The print method is called for its side effect, i.e. printing the summary.

    +
    +
    +

    References

    +

    FOCUS (2006) “Guidance Document on Estimating Persistence and Degradation Kinetics from Environmental Fate Studies on Pesticides in -EU Registration” Report of the FOCUS Work Group on Degradation Kinetics, +EU Registration” Report of the FOCUS Work Group on Degradation Kinetics, EC Document Reference Sanco/10058/2005 version 2.0, 434 pp, -http://esdac.jrc.ec.europa.eu/projects/degradation-kinetics

    -

    Author

    - +http://esdac.jrc.ec.europa.eu/projects/degradation-kinetics

    +
    +
    +

    Author

    Johannes Ranke

    +
    -

    Examples

    - 
    
-  summary(mkinfit(mkinmod(parent = mkinsub("SFO")), FOCUS_2006_A, quiet = TRUE))
-
    #> mkin version used for fitting:    1.0.3.9000 
-#> R version used for fitting:       4.0.3 
-#> Date of fit:     Mon Feb 15 17:13:09 2021 
-#> Date of summary: Mon Feb 15 17:13:09 2021 
-#> 
-#> Equations:
-#> d_parent/dt = - k_parent * parent
-#> 
-#> Model predictions using solution type analytical 
-#> 
-#> Fitted using 131 model solutions performed in 0.027 s
-#> 
-#> Error model: Constant variance 
-#> 
-#> Error model algorithm: OLS 
-#> 
-#> Starting values for parameters to be optimised:
-#>           value   type
-#> parent_0 101.24  state
-#> k_parent   0.10 deparm
-#> 
-#> Starting values for the transformed parameters actually optimised:
-#>                   value lower upper
-#> parent_0     101.240000  -Inf   Inf
-#> log_k_parent  -2.302585  -Inf   Inf
-#> 
-#> Fixed parameter values:
-#> None
-#> 
-#> Results:
-#> 
-#>        AIC     BIC    logLik
-#>   55.28197 55.5203 -24.64099
-#> 
-#> Optimised, transformed parameters with symmetric confidence intervals:
-#>              Estimate Std. Error  Lower   Upper
-#> parent_0      109.200    3.70400 99.630 118.700
-#> log_k_parent   -3.291    0.09176 -3.527  -3.055
-#> sigma           5.266    1.31600  1.882   8.649
-#> 
-#> Parameter correlation:
-#>               parent_0 log_k_parent     sigma
-#> parent_0     1.000e+00    5.428e-01 1.642e-07
-#> log_k_parent 5.428e-01    1.000e+00 2.507e-07
-#> sigma        1.642e-07    2.507e-07 1.000e+00
-#> 
-#> Backtransformed parameters:
-#> Confidence intervals for internally transformed parameters are asymmetric.
-#> t-test (unrealistically) based on the assumption of normal distribution
-#> for estimators of untransformed parameters.
-#>           Estimate t value    Pr(>t)   Lower     Upper
-#> parent_0 109.20000   29.47 4.218e-07 99.6300 118.70000
-#> k_parent   0.03722   10.90 5.650e-05  0.0294   0.04712
-#> sigma      5.26600    4.00 5.162e-03  1.8820   8.64900
-#> 
-#> FOCUS Chi2 error levels in percent:
-#>          err.min n.optim df
-#> All data   8.385       2  6
-#> parent     8.385       2  6
-#> 
-#> Estimated disappearance times:
-#>         DT50  DT90
-#> parent 18.62 61.87
-#> 
-#> Data:
-#>  time variable observed predicted residual
-#>     0   parent   101.24   109.153  -7.9132
-#>     3   parent    99.27    97.622   1.6484
-#>     7   parent    90.11    84.119   5.9913
-#>    14   parent    72.19    64.826   7.3641
-#>    30   parent    29.71    35.738  -6.0283
-#>    62   parent     5.98    10.862  -4.8818
-#>    90   parent     1.54     3.831  -2.2911
-#>   118   parent     0.39     1.351  -0.9613
    
-
    +
    +

    Examples

    + 
    
+  summary(mkinfit("SFO", FOCUS_2006_A, quiet = TRUE))
+#> mkin version used for fitting:    1.1.2 
+#> R version used for fitting:       4.2.1 
+#> Date of fit:     Wed Aug 10 15:27:28 2022 
+#> Date of summary: Wed Aug 10 15:27:28 2022 
+#> 
+#> Equations:
+#> d_parent/dt = - k_parent * parent
+#> 
+#> Model predictions using solution type analytical 
+#> 
+#> Fitted using 131 model solutions performed in 0.029 s
+#> 
+#> Error model: Constant variance 
+#> 
+#> Error model algorithm: OLS 
+#> 
+#> Starting values for parameters to be optimised:
+#>           value   type
+#> parent_0 101.24  state
+#> k_parent   0.10 deparm
+#> 
+#> Starting values for the transformed parameters actually optimised:
+#>                   value lower upper
+#> parent_0     101.240000  -Inf   Inf
+#> log_k_parent  -2.302585  -Inf   Inf
+#> 
+#> Fixed parameter values:
+#> None
+#> 
+#> Results:
+#> 
+#>        AIC     BIC    logLik
+#>   55.28197 55.5203 -24.64099
+#> 
+#> Optimised, transformed parameters with symmetric confidence intervals:
+#>              Estimate Std. Error  Lower   Upper
+#> parent_0      109.200    3.70400 99.630 118.700
+#> log_k_parent   -3.291    0.09176 -3.527  -3.055
+#> sigma           5.266    1.31600  1.882   8.649
+#> 
+#> Parameter correlation:
+#>               parent_0 log_k_parent     sigma
+#> parent_0     1.000e+00    5.428e-01 1.642e-07
+#> log_k_parent 5.428e-01    1.000e+00 2.507e-07
+#> sigma        1.642e-07    2.507e-07 1.000e+00
+#> 
+#> Backtransformed parameters:
+#> Confidence intervals for internally transformed parameters are asymmetric.
+#> t-test (unrealistically) based on the assumption of normal distribution
+#> for estimators of untransformed parameters.
+#>           Estimate t value    Pr(>t)   Lower     Upper
+#> parent_0 109.20000   29.47 4.218e-07 99.6300 118.70000
+#> k_parent   0.03722   10.90 5.650e-05  0.0294   0.04712
+#> sigma      5.26600    4.00 5.162e-03  1.8820   8.64900
+#> 
+#> FOCUS Chi2 error levels in percent:
+#>          err.min n.optim df
+#> All data   8.385       2  6
+#> parent     8.385       2  6
+#> 
+#> Estimated disappearance times:
+#>         DT50  DT90
+#> parent 18.62 61.87
+#> 
+#> Data:
+#>  time variable observed predicted residual
+#>     0   parent   101.24   109.153  -7.9132
+#>     3   parent    99.27    97.622   1.6484
+#>     7   parent    90.11    84.119   5.9913
+#>    14   parent    72.19    64.826   7.3641
+#>    30   parent    29.71    35.738  -6.0283
+#>    62   parent     5.98    10.862  -4.8818
+#>    90   parent     1.54     3.831  -2.2911
+#>   118   parent     0.39     1.351  -0.9613
+
+
    +
    +
    -
    -
    -

    Developed by Johannes Ranke.

    +
    +

    Developed by Johannes Ranke.

    -

    Site built with pkgdown 1.6.1.

    +

    Site built with pkgdown 2.0.6.

    -
    -
    +
    - - + + diff --git a/docs/dev/reference/summary.mmkin.html b/docs/dev/reference/summary.mmkin.html new file mode 100644 index 00000000..e916ba24 --- /dev/null +++ b/docs/dev/reference/summary.mmkin.html @@ -0,0 +1,173 @@ + +Summary method for class "mmkin" — summary.mmkin • mkin + + +
    +
    + + + +
    +
    +
    +

    Summary method for class "mmkin"

    + Source: R/summary.mmkin.R + +
    + +
    +

    Shows convergence information on the mkinfit objects contained in the object +and gives an overview of ill-defined parameters calculated by illparms.

    +
    + +
    + 
    # S3 method for mmkin
+summary(object, conf.level = 0.95, ...)
+
+# S3 method for summary.mmkin
+print(x, digits = max(3, getOption("digits") - 3), ...)
    +
    + +
    +

    Arguments

    +
    object
    +

    an object of class mmkin

    + + +
    conf.level
    +

    confidence level for testing parameters

    + + +
    ...
    +

    optional arguments passed to methods like print.

    + + +
    x
    +

    an object of class summary.mmkin.

    + + +
    digits
    +

    number of digits to use for printing

    + +
    + +
    +

    Examples

    + 
    
+fits <- mmkin(
+  c("SFO", "FOMC"),
+  list("FOCUS A" = FOCUS_2006_A,
+       "FOCUS C" = FOCUS_2006_C),
+  quiet = TRUE, cores = 1)
+  summary(fits)
+#> Error model: Constant variance 
+#> Fitted in 0.802 s
+#> 
+#> Convergence:
+#>       dataset
+#> model  FOCUS A FOCUS C
+#>   SFO  OK      OK     
+#>   FOMC OK      OK     
+#> 
+#> OK: No warnings
+#> 
+#> Ill-defined parameters:
+#>       dataset
+#> model  FOCUS A                      FOCUS C
+#>   SFO                                      
+#>   FOMC parent_0, alpha, beta, sigma        
+
+
    +
    +
    + +
    + + +
    +

    Developed by Johannes Ranke.

    +
    + +
    +

    Site built with pkgdown 2.0.6.

    +
    + +
    + + + + + + + + diff --git a/docs/dev/reference/summary.nlme.mmkin.html b/docs/dev/reference/summary.nlme.mmkin.html index b2f6624a..067efcfe 100644 --- a/docs/dev/reference/summary.nlme.mmkin.html +++ b/docs/dev/reference/summary.nlme.mmkin.html @@ -1,71 +1,16 @@ - - - - - - - -Summary method for class "nlme.mmkin" — summary.nlme.mmkin • mkin - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -Summary method for class "nlme.mmkin" — summary.nlme.mmkin • mkin - - - - - - - - - - + + - - - -
    -
    - -
    - -
    +

    Summary method for class "nlme.mmkin"

    - Source: R/summary.nlme.mmkin.R + Source: R/summary.nlme.mmkin.R
    @@ -155,295 +94,320 @@ endpoints such as formation fractions and DT50 values. Optionally (default is FALSE), the data are listed in full.

    - 
    # S3 method for nlme.mmkin
-summary(
-  object,
-  data = FALSE,
-  verbose = FALSE,
-  distimes = TRUE,
-  alpha = 0.05,
-  ...
-)
-
-# S3 method for summary.nlme.mmkin
-print(x, digits = max(3, getOption("digits") - 3), verbose = x$verbose, ...)
    - -

    Arguments

    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    object

    an object of class nlme.mmkin

    data

    logical, indicating whether the full data should be included in -the summary.

    verbose

    Should the summary be verbose?

    distimes

    logical, indicating whether DT50 and DT90 values should be -included.

    alpha

    error level for confidence interval estimation from the t -distribution

    ...

    optional arguments passed to methods like print.

    x

    an object of class summary.nlme.mmkin

    digits

    Number of digits to use for printing

    - -

    Value

    - -

    The summary function returns a list based on the nlme object +

    + 
    # S3 method for nlme.mmkin
+summary(
+  object,
+  data = FALSE,
+  verbose = FALSE,
+  distimes = TRUE,
+  alpha = 0.05,
+  ...
+)
+
+# S3 method for summary.nlme.mmkin
+print(x, digits = max(3, getOption("digits") - 3), verbose = x$verbose, ...)
    +
    + +
    +

    Arguments

    +
    object
    +

    an object of class nlme.mmkin

    + + +
    data
    +

    logical, indicating whether the full data should be included in +the summary.

    + + +
    verbose
    +

    Should the summary be verbose?

    + + +
    distimes
    +

    logical, indicating whether DT50 and DT90 values should be +included.

    + + +
    alpha
    +

    error level for confidence interval estimation from the t +distribution

    + + +
    ...
    +

    optional arguments passed to methods like print.

    + + +
    x
    +

    an object of class summary.nlme.mmkin

    + + +
    digits
    +

    Number of digits to use for printing

    + +
    +
    +

    Value

    + + +

    The summary function returns a list based on the nlme object obtained in the fit, with at least the following additional components

    -
    nlmeversion, mkinversion, Rversion

    The nlme, mkin and R versions used

    -
    date.fit, date.summary

    The dates where the fit and the summary were +

    nlmeversion, mkinversion, Rversion
    +

    The nlme, mkin and R versions used

    + +
    date.fit, date.summary
    +

    The dates where the fit and the summary were produced

    -
    diffs

    The differential equations used in the degradation model

    -
    use_of_ff

    Was maximum or minimum use made of formation fractions

    -
    data

    The data

    -
    confint_trans

    Transformed parameters as used in the optimisation, with confidence intervals

    -
    confint_back

    Backtransformed parameters, with confidence intervals if available

    -
    ff

    The estimated formation fractions derived from the fitted + +

    diffs
    +

    The differential equations used in the degradation model

    + +
    use_of_ff
    +

    Was maximum or minimum use made of formation fractions

    + +
    data
    +

    The data

    + +
    confint_trans
    +

    Transformed parameters as used in the optimisation, with confidence intervals

    + +
    confint_back
    +

    Backtransformed parameters, with confidence intervals if available

    + +
    ff
    +

    The estimated formation fractions derived from the fitted model.

    -
    distimes

    The DT50 and DT90 values for each observed variable.

    -
    SFORB

    If applicable, eigenvalues of SFORB components of the model.

    -The print method is called for its side effect, i.e. printing the summary. -

    Author

    +
    distimes
    +

    The DT50 and DT90 values for each observed variable.

    + +
    SFORB
    +

    If applicable, eigenvalues of SFORB components of the model.

    +

    The print method is called for its side effect, i.e. printing the summary.

    +
    +
    +

    Author

    Johannes Ranke for the mkin specific parts José Pinheiro and Douglas Bates for the components inherited from nlme

    +
    -

    Examples

    - 
    
-# Generate five datasets following SFO kinetics
-sampling_times = c(0, 1, 3, 7, 14, 28, 60, 90, 120)
-dt50_sfo_in_pop <- 50
-k_in_pop <- log(2) / dt50_sfo_in_pop
-set.seed(1234)
-k_in <- rlnorm(5, log(k_in_pop), 0.5)
-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)
-names(ds_sfo_mean) <- paste("ds", 1:5)
-
-set.seed(12345)
-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]]
-})
-
-# Evaluate using mmkin and nlme
-library(nlme)
-f_mmkin <- mmkin("SFO", ds_sfo_syn, quiet = TRUE, error_model = "tc", cores = 1)
-
    #> Warning: Optimisation did not converge:
-#> iteration limit reached without convergence (10)
    f_nlme <- nlme(f_mmkin)
-
    #> Warning: Iteration 4, LME step: nlminb() did not converge (code = 1). PORT message: false convergence (8)
    summary(f_nlme, data = TRUE)
-
    #> nlme version used for fitting:      3.1.152 
-#> mkin version used for pre-fitting:  1.0.3.9000 
-#> R version used for fitting:         4.0.3 
-#> Date of fit:     Mon Feb 15 17:13:12 2021 
-#> Date of summary: Mon Feb 15 17:13:12 2021 
-#> 
-#> Equations:
-#> d_parent/dt = - k_parent * parent
-#> 
-#> Data:
-#> 90 observations of 1 variable(s) grouped in 5 datasets
-#> 
-#> Model predictions using solution type analytical 
-#> 
-#> Fitted in 0.544 s using 4 iterations
-#> 
-#> Variance model: Two-component variance function 
-#> 
-#> Mean of starting values for individual parameters:
-#>     parent_0 log_k_parent 
-#>      101.569       -4.454 
-#> 
-#> Fixed degradation parameter values:
-#> None
-#> 
-#> Results:
-#> 
-#>     AIC   BIC logLik
-#>   584.5 599.5 -286.2
-#> 
-#> Optimised, transformed parameters with symmetric confidence intervals:
-#>               lower    est.   upper
-#> parent_0     99.371 101.592 103.814
-#> log_k_parent -4.973  -4.449  -3.926
-#> 
-#> Correlation: 
-#>              prnt_0
-#> log_k_parent 0.051 
-#> 
-#> Random effects:
-#>  Formula: list(parent_0 ~ 1, log_k_parent ~ 1)
-#>  Level: ds
-#>  Structure: Diagonal
-#>          parent_0 log_k_parent Residual
-#> StdDev: 6.924e-05       0.5863        1
-#> 
-#> Variance function:
-#>  Structure: Constant plus proportion of variance covariate
-#>  Formula: ~fitted(.) 
-#>  Parameter estimates:
-#>        const         prop 
-#> 0.0001208853 0.0789968036 
-#> 
-#> Backtransformed parameters with asymmetric confidence intervals:
-#>              lower      est.     upper
-#> parent_0 99.370882 101.59243 103.81398
-#> k_parent  0.006923   0.01168   0.01972
-#> 
-#> Estimated disappearance times:
-#>         DT50  DT90
-#> parent 59.32 197.1
-#> 
-#> Data:
-#>    ds   name time observed predicted  residual    std standardized
-#>  ds 1 parent    0    104.1   101.592   2.50757 8.0255     0.312451
-#>  ds 1 parent    0    105.0   101.592   3.40757 8.0255     0.424594
-#>  ds 1 parent    1     98.5   100.796  -2.29571 7.9625    -0.288313
-#>  ds 1 parent    1     96.1   100.796  -4.69571 7.9625    -0.589725
-#>  ds 1 parent    3    101.9    99.221   2.67904 7.8381     0.341796
-#>  ds 1 parent    3     85.2    99.221 -14.02096 7.8381    -1.788812
-#>  ds 1 parent    7     99.1    96.145   2.95512 7.5951     0.389081
-#>  ds 1 parent    7     93.0    96.145  -3.14488 7.5951    -0.414065
-#>  ds 1 parent   14     88.1    90.989  -2.88944 7.1879    -0.401987
-#>  ds 1 parent   14     84.1    90.989  -6.88944 7.1879    -0.958480
-#>  ds 1 parent   28     80.2    81.493  -1.29305 6.4377    -0.200857
-#>  ds 1 parent   28     91.3    81.493   9.80695 6.4377     1.523364
-#>  ds 1 parent   60     65.1    63.344   1.75642 5.0039     0.351008
-#>  ds 1 parent   60     65.8    63.344   2.45642 5.0039     0.490898
-#>  ds 1 parent   90     47.8    50.018  -2.21764 3.9512    -0.561252
-#>  ds 1 parent   90     53.5    50.018   3.48236 3.9512     0.881335
-#>  ds 1 parent  120     37.6    39.495  -1.89515 3.1200    -0.607423
-#>  ds 1 parent  120     39.3    39.495  -0.19515 3.1200    -0.062549
-#>  ds 2 parent    0    107.9   101.592   6.30757 8.0255     0.785943
-#>  ds 2 parent    0    102.1   101.592   0.50757 8.0255     0.063245
-#>  ds 2 parent    1    103.8   100.058   3.74159 7.9043     0.473361
-#>  ds 2 parent    1    108.6   100.058   8.54159 7.9043     1.080626
-#>  ds 2 parent    3     91.0    97.060  -6.05952 7.6674    -0.790297
-#>  ds 2 parent    3     84.9    97.060 -12.15952 7.6674    -1.585874
-#>  ds 2 parent    7     79.3    91.329 -12.02867 7.2147    -1.667251
-#>  ds 2 parent    7    100.9    91.329   9.57133 7.2147     1.326647
-#>  ds 2 parent   14     77.3    82.102  -4.80185 6.4858    -0.740366
-#>  ds 2 parent   14     83.5    82.102   1.39815 6.4858     0.215571
-#>  ds 2 parent   28     66.8    66.351   0.44945 5.2415     0.085748
-#>  ds 2 parent   28     63.3    66.351  -3.05055 5.2415    -0.582002
-#>  ds 2 parent   60     40.8    40.775   0.02474 3.2211     0.007679
-#>  ds 2 parent   60     44.8    40.775   4.02474 3.2211     1.249485
-#>  ds 2 parent   90     27.8    25.832   1.96762 2.0407     0.964198
-#>  ds 2 parent   90     27.0    25.832   1.16762 2.0407     0.572171
-#>  ds 2 parent  120     15.2    16.366  -1.16561 1.2928    -0.901595
-#>  ds 2 parent  120     15.5    16.366  -0.86561 1.2928    -0.669547
-#>  ds 3 parent    0     97.7   101.592  -3.89243 8.0255    -0.485009
-#>  ds 3 parent    0     88.2   101.592 -13.39243 8.0255    -1.668739
-#>  ds 3 parent    1    109.9    99.218  10.68196 7.8379     1.362858
-#>  ds 3 parent    1     97.8    99.218  -1.41804 7.8379    -0.180921
-#>  ds 3 parent    3    100.5    94.634   5.86555 7.4758     0.784603
-#>  ds 3 parent    3     77.4    94.634 -17.23445 7.4758    -2.305360
-#>  ds 3 parent    7     78.3    86.093  -7.79273 6.8011    -1.145813
-#>  ds 3 parent    7     90.3    86.093   4.20727 6.8011     0.618620
-#>  ds 3 parent   14     76.0    72.958   3.04222 5.7634     0.527848
-#>  ds 3 parent   14     79.1    72.958   6.14222 5.7634     1.065722
-#>  ds 3 parent   28     46.0    52.394  -6.39404 4.1390    -1.544842
-#>  ds 3 parent   28     53.4    52.394   1.00596 4.1390     0.243046
-#>  ds 3 parent   60     25.1    24.582   0.51786 1.9419     0.266676
-#>  ds 3 parent   60     21.4    24.582  -3.18214 1.9419    -1.638664
-#>  ds 3 parent   90     11.0    12.092  -1.09202 0.9552    -1.143199
-#>  ds 3 parent   90     14.2    12.092   2.10798 0.9552     2.206776
-#>  ds 3 parent  120      5.8     5.948  -0.14810 0.4699    -0.315178
-#>  ds 3 parent  120      6.1     5.948   0.15190 0.4699     0.323282
-#>  ds 4 parent    0     95.3   101.592  -6.29243 8.0255    -0.784057
-#>  ds 4 parent    0    102.0   101.592   0.40757 8.0255     0.050784
-#>  ds 4 parent    1    104.4   101.125   3.27549 7.9885     0.410025
-#>  ds 4 parent    1    105.4   101.125   4.27549 7.9885     0.535205
-#>  ds 4 parent    3    113.7   100.195  13.50487 7.9151     1.706218
-#>  ds 4 parent    3     82.3   100.195 -17.89513 7.9151    -2.260886
-#>  ds 4 parent    7     98.1    98.362  -0.26190 7.7703    -0.033706
-#>  ds 4 parent    7     87.8    98.362 -10.56190 7.7703    -1.359270
-#>  ds 4 parent   14     97.9    95.234   2.66590 7.5232     0.354357
-#>  ds 4 parent   14    104.8    95.234   9.56590 7.5232     1.271521
-#>  ds 4 parent   28     85.0    89.274  -4.27372 7.0523    -0.606001
-#>  ds 4 parent   28     77.2    89.274 -12.07372 7.0523    -1.712017
-#>  ds 4 parent   60     82.2    77.013   5.18661 6.0838     0.852526
-#>  ds 4 parent   60     86.1    77.013   9.08661 6.0838     1.493571
-#>  ds 4 parent   90     70.5    67.053   3.44692 5.2970     0.650733
-#>  ds 4 parent   90     61.7    67.053  -5.35308 5.2970    -1.010591
-#>  ds 4 parent  120     60.0    58.381   1.61905 4.6119     0.351058
-#>  ds 4 parent  120     56.4    58.381  -1.98095 4.6119    -0.429530
-#>  ds 5 parent    0     92.6   101.592  -8.99243 8.0255    -1.120485
-#>  ds 5 parent    0    116.5   101.592  14.90757 8.0255     1.857531
-#>  ds 5 parent    1    108.0    99.914   8.08560 7.8929     1.024413
-#>  ds 5 parent    1    104.9    99.914   4.98560 7.8929     0.631655
-#>  ds 5 parent    3    100.5    96.641   3.85898 7.6343     0.505477
-#>  ds 5 parent    3     89.5    96.641  -7.14102 7.6343    -0.935382
-#>  ds 5 parent    7     91.7    90.412   1.28752 7.1423     0.180267
-#>  ds 5 parent    7     95.1    90.412   4.68752 7.1423     0.656304
-#>  ds 5 parent   14     82.2    80.463   1.73715 6.3563     0.273295
-#>  ds 5 parent   14     84.5    80.463   4.03715 6.3563     0.635141
-#>  ds 5 parent   28     60.5    63.728  -3.22788 5.0343    -0.641178
-#>  ds 5 parent   28     72.8    63.728   9.07212 5.0343     1.802062
-#>  ds 5 parent   60     38.3    37.399   0.90061 2.9544     0.304835
-#>  ds 5 parent   60     40.7    37.399   3.30061 2.9544     1.117174
-#>  ds 5 parent   90     22.5    22.692  -0.19165 1.7926    -0.106913
-#>  ds 5 parent   90     20.8    22.692  -1.89165 1.7926    -1.055273
-#>  ds 5 parent  120     13.4    13.768  -0.36790 1.0876    -0.338259
-#>  ds 5 parent  120     13.8    13.768   0.03210 1.0876     0.029517
    
-
    +
    +

    Examples

    + 
    
+# Generate five datasets following SFO kinetics
+sampling_times = c(0, 1, 3, 7, 14, 28, 60, 90, 120)
+dt50_sfo_in_pop <- 50
+k_in_pop <- log(2) / dt50_sfo_in_pop
+set.seed(1234)
+k_in <- rlnorm(5, log(k_in_pop), 0.5)
+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)
+names(ds_sfo_mean) <- paste("ds", 1:5)
+
+set.seed(12345)
+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]]
+})
+
+# \dontrun{
+# Evaluate using mmkin and nlme
+library(nlme)
+f_mmkin <- mmkin("SFO", ds_sfo_syn, quiet = TRUE, error_model = "tc", cores = 1)
+#> Warning: Optimisation did not converge:
+#> iteration limit reached without convergence (10)
+f_nlme <- nlme(f_mmkin)
+#> Warning: Iteration 4, LME step: nlminb() did not converge (code = 1). PORT message: false convergence (8)
+summary(f_nlme, data = TRUE)
+#> nlme version used for fitting:      3.1.158 
+#> mkin version used for pre-fitting:  1.1.2 
+#> R version used for fitting:         4.2.1 
+#> Date of fit:     Wed Aug 10 15:27:32 2022 
+#> Date of summary: Wed Aug 10 15:27:32 2022 
+#> 
+#> Equations:
+#> d_parent/dt = - k_parent * parent
+#> 
+#> Data:
+#> 90 observations of 1 variable(s) grouped in 5 datasets
+#> 
+#> Model predictions using solution type analytical 
+#> 
+#> Fitted in 0.534 s using 4 iterations
+#> 
+#> Variance model: Two-component variance function 
+#> 
+#> Mean of starting values for individual parameters:
+#>     parent_0 log_k_parent 
+#>      101.569       -4.454 
+#> 
+#> Fixed degradation parameter values:
+#> None
+#> 
+#> Results:
+#> 
+#>     AIC   BIC logLik
+#>   584.5 599.5 -286.2
+#> 
+#> Optimised, transformed parameters with symmetric confidence intervals:
+#>               lower    est.   upper
+#> parent_0     99.371 101.592 103.814
+#> log_k_parent -4.973  -4.449  -3.926
+#> 
+#> Correlation: 
+#>              parnt_0
+#> log_k_parent 0.0507 
+#> 
+#> Random effects:
+#>  Formula: list(parent_0 ~ 1, log_k_parent ~ 1)
+#>  Level: ds
+#>  Structure: Diagonal
+#>          parent_0 log_k_parent Residual
+#> StdDev: 6.924e-05       0.5863        1
+#> 
+#> Variance function:
+#>  Structure: Constant plus proportion of variance covariate
+#>  Formula: ~fitted(.) 
+#>  Parameter estimates:
+#>        const         prop 
+#> 0.0001208853 0.0789968036 
+#> 
+#> Backtransformed parameters with asymmetric confidence intervals:
+#>              lower      est.     upper
+#> parent_0 99.370882 101.59243 103.81398
+#> k_parent  0.006923   0.01168   0.01972
+#> 
+#> Estimated disappearance times:
+#>         DT50  DT90
+#> parent 59.32 197.1
+#> 
+#> Data:
+#>    ds   name time observed predicted  residual    std standardized
+#>  ds 1 parent    0    104.1   101.592   2.50757 8.0255     0.312451
+#>  ds 1 parent    0    105.0   101.592   3.40757 8.0255     0.424594
+#>  ds 1 parent    1     98.5   100.796  -2.29571 7.9625    -0.288313
+#>  ds 1 parent    1     96.1   100.796  -4.69571 7.9625    -0.589725
+#>  ds 1 parent    3    101.9    99.221   2.67904 7.8381     0.341796
+#>  ds 1 parent    3     85.2    99.221 -14.02096 7.8381    -1.788812
+#>  ds 1 parent    7     99.1    96.145   2.95512 7.5951     0.389081
+#>  ds 1 parent    7     93.0    96.145  -3.14488 7.5951    -0.414065
+#>  ds 1 parent   14     88.1    90.989  -2.88944 7.1879    -0.401987
+#>  ds 1 parent   14     84.1    90.989  -6.88944 7.1879    -0.958480
+#>  ds 1 parent   28     80.2    81.493  -1.29305 6.4377    -0.200857
+#>  ds 1 parent   28     91.3    81.493   9.80695 6.4377     1.523364
+#>  ds 1 parent   60     65.1    63.344   1.75642 5.0039     0.351008
+#>  ds 1 parent   60     65.8    63.344   2.45642 5.0039     0.490898
+#>  ds 1 parent   90     47.8    50.018  -2.21764 3.9512    -0.561252
+#>  ds 1 parent   90     53.5    50.018   3.48236 3.9512     0.881335
+#>  ds 1 parent  120     37.6    39.495  -1.89515 3.1200    -0.607423
+#>  ds 1 parent  120     39.3    39.495  -0.19515 3.1200    -0.062549
+#>  ds 2 parent    0    107.9   101.592   6.30757 8.0255     0.785943
+#>  ds 2 parent    0    102.1   101.592   0.50757 8.0255     0.063245
+#>  ds 2 parent    1    103.8   100.058   3.74159 7.9043     0.473361
+#>  ds 2 parent    1    108.6   100.058   8.54159 7.9043     1.080626
+#>  ds 2 parent    3     91.0    97.060  -6.05952 7.6674    -0.790297
+#>  ds 2 parent    3     84.9    97.060 -12.15952 7.6674    -1.585874
+#>  ds 2 parent    7     79.3    91.329 -12.02867 7.2147    -1.667251
+#>  ds 2 parent    7    100.9    91.329   9.57133 7.2147     1.326647
+#>  ds 2 parent   14     77.3    82.102  -4.80185 6.4858    -0.740366
+#>  ds 2 parent   14     83.5    82.102   1.39815 6.4858     0.215571
+#>  ds 2 parent   28     66.8    66.351   0.44945 5.2415     0.085748
+#>  ds 2 parent   28     63.3    66.351  -3.05055 5.2415    -0.582002
+#>  ds 2 parent   60     40.8    40.775   0.02474 3.2211     0.007679
+#>  ds 2 parent   60     44.8    40.775   4.02474 3.2211     1.249485
+#>  ds 2 parent   90     27.8    25.832   1.96762 2.0407     0.964198
+#>  ds 2 parent   90     27.0    25.832   1.16762 2.0407     0.572171
+#>  ds 2 parent  120     15.2    16.366  -1.16561 1.2928    -0.901595
+#>  ds 2 parent  120     15.5    16.366  -0.86561 1.2928    -0.669547
+#>  ds 3 parent    0     97.7   101.592  -3.89243 8.0255    -0.485009
+#>  ds 3 parent    0     88.2   101.592 -13.39243 8.0255    -1.668739
+#>  ds 3 parent    1    109.9    99.218  10.68196 7.8379     1.362858
+#>  ds 3 parent    1     97.8    99.218  -1.41804 7.8379    -0.180921
+#>  ds 3 parent    3    100.5    94.634   5.86555 7.4758     0.784603
+#>  ds 3 parent    3     77.4    94.634 -17.23445 7.4758    -2.305360
+#>  ds 3 parent    7     78.3    86.093  -7.79273 6.8011    -1.145813
+#>  ds 3 parent    7     90.3    86.093   4.20727 6.8011     0.618620
+#>  ds 3 parent   14     76.0    72.958   3.04222 5.7634     0.527848
+#>  ds 3 parent   14     79.1    72.958   6.14222 5.7634     1.065722
+#>  ds 3 parent   28     46.0    52.394  -6.39404 4.1390    -1.544842
+#>  ds 3 parent   28     53.4    52.394   1.00596 4.1390     0.243046
+#>  ds 3 parent   60     25.1    24.582   0.51786 1.9419     0.266676
+#>  ds 3 parent   60     21.4    24.582  -3.18214 1.9419    -1.638664
+#>  ds 3 parent   90     11.0    12.092  -1.09202 0.9552    -1.143199
+#>  ds 3 parent   90     14.2    12.092   2.10798 0.9552     2.206776
+#>  ds 3 parent  120      5.8     5.948  -0.14810 0.4699    -0.315178
+#>  ds 3 parent  120      6.1     5.948   0.15190 0.4699     0.323282
+#>  ds 4 parent    0     95.3   101.592  -6.29243 8.0255    -0.784057
+#>  ds 4 parent    0    102.0   101.592   0.40757 8.0255     0.050784
+#>  ds 4 parent    1    104.4   101.125   3.27549 7.9885     0.410025
+#>  ds 4 parent    1    105.4   101.125   4.27549 7.9885     0.535205
+#>  ds 4 parent    3    113.7   100.195  13.50487 7.9151     1.706218
+#>  ds 4 parent    3     82.3   100.195 -17.89513 7.9151    -2.260886
+#>  ds 4 parent    7     98.1    98.362  -0.26190 7.7703    -0.033706
+#>  ds 4 parent    7     87.8    98.362 -10.56190 7.7703    -1.359270
+#>  ds 4 parent   14     97.9    95.234   2.66590 7.5232     0.354357
+#>  ds 4 parent   14    104.8    95.234   9.56590 7.5232     1.271521
+#>  ds 4 parent   28     85.0    89.274  -4.27372 7.0523    -0.606001
+#>  ds 4 parent   28     77.2    89.274 -12.07372 7.0523    -1.712017
+#>  ds 4 parent   60     82.2    77.013   5.18661 6.0838     0.852526
+#>  ds 4 parent   60     86.1    77.013   9.08661 6.0838     1.493571
+#>  ds 4 parent   90     70.5    67.053   3.44692 5.2970     0.650733
+#>  ds 4 parent   90     61.7    67.053  -5.35308 5.2970    -1.010591
+#>  ds 4 parent  120     60.0    58.381   1.61905 4.6119     0.351058
+#>  ds 4 parent  120     56.4    58.381  -1.98095 4.6119    -0.429530
+#>  ds 5 parent    0     92.6   101.592  -8.99243 8.0255    -1.120485
+#>  ds 5 parent    0    116.5   101.592  14.90757 8.0255     1.857531
+#>  ds 5 parent    1    108.0    99.914   8.08560 7.8929     1.024413
+#>  ds 5 parent    1    104.9    99.914   4.98560 7.8929     0.631655
+#>  ds 5 parent    3    100.5    96.641   3.85898 7.6343     0.505477
+#>  ds 5 parent    3     89.5    96.641  -7.14102 7.6343    -0.935382
+#>  ds 5 parent    7     91.7    90.412   1.28752 7.1423     0.180267
+#>  ds 5 parent    7     95.1    90.412   4.68752 7.1423     0.656304
+#>  ds 5 parent   14     82.2    80.463   1.73715 6.3563     0.273295
+#>  ds 5 parent   14     84.5    80.463   4.03715 6.3563     0.635141
+#>  ds 5 parent   28     60.5    63.728  -3.22788 5.0343    -0.641178
+#>  ds 5 parent   28     72.8    63.728   9.07212 5.0343     1.802062
+#>  ds 5 parent   60     38.3    37.399   0.90061 2.9544     0.304835
+#>  ds 5 parent   60     40.7    37.399   3.30061 2.9544     1.117174
+#>  ds 5 parent   90     22.5    22.692  -0.19165 1.7926    -0.106913
+#>  ds 5 parent   90     20.8    22.692  -1.89165 1.7926    -1.055273
+#>  ds 5 parent  120     13.4    13.768  -0.36790 1.0876    -0.338259
+#>  ds 5 parent  120     13.8    13.768   0.03210 1.0876     0.029517
+# }
+
+
    +
    +
    -
    -
    -

    Developed by Johannes Ranke.

    +
    +

    Developed by Johannes Ranke.

    -

    Site built with pkgdown 1.6.1.

    +

    Site built with pkgdown 2.0.6.

    -
    -
    +
    - - + + diff --git a/docs/dev/reference/summary.saem.mmkin.html b/docs/dev/reference/summary.saem.mmkin.html index febf25b4..6cab9c5f 100644 --- a/docs/dev/reference/summary.saem.mmkin.html +++ b/docs/dev/reference/summary.saem.mmkin.html @@ -21,7 +21,7 @@ endpoints such as formation fractions and DT50 values. Optionally mkin - 1.1.0 + 1.1.2 @@ -30,7 +30,7 @@ endpoints such as formation fractions and DT50 values. Optionally Functions and data
  • -
  • Example evaluation of FOCUS Laboratory Data L1 to L3
    • +
  • + Example evaluations of dimethenamid data from 2018 with nonlinear mixed-effects models +
  • Example evaluation of FOCUS Example Dataset Z
    • @@ -92,59 +95,88 @@ endpoints such as formation fractions and DT50 values. Optionally
    - 
    # S3 method for saem.mmkin
-summary(object, data = FALSE, verbose = FALSE, distimes = TRUE, ...)
-
-# S3 method for summary.saem.mmkin
-print(x, digits = max(3, getOption("digits") - 3), verbose = x$verbose, ...)
    + 
    # S3 method for saem.mmkin
+summary(object, data = FALSE, verbose = FALSE, distimes = TRUE, ...)
+
+# S3 method for summary.saem.mmkin
+print(x, digits = max(3, getOption("digits") - 3), verbose = x$verbose, ...)

    Arguments

    object

    an object of class saem.mmkin

    + +
    data

    logical, indicating whether the full data should be included in the summary.

    + +
    verbose

    Should the summary be verbose?

    + +
    distimes

    logical, indicating whether DT50 and DT90 values should be included.

    + +
    ...

    optional arguments passed to methods like print.

    + +
    x

    an object of class summary.saem.mmkin

    + +
    digits

    Number of digits to use for printing

    +

    Value

    -

    The summary function returns a list based on the saemix::SaemixObjectobtained in the fit, with at least the following additional components

    + + +

    The summary function returns a list based on the saemix::SaemixObject

    + + +

    obtained in the fit, with at least the following additional components

    saemixversion, mkinversion, Rversion

    The saemix, mkin and R versions used

    +
    date.fit, date.summary

    The dates where the fit and the summary were produced

    +
    diffs

    The differential equations used in the degradation model

    +
    use_of_ff

    Was maximum or minimum use made of formation fractions

    +
    data

    The data

    +
    confint_trans

    Transformed parameters as used in the optimisation, with confidence intervals

    +
    confint_back

    Backtransformed parameters, with confidence intervals if available

    +
    confint_errmod

    Error model parameters with confidence intervals

    +
    ff

    The estimated formation fractions derived from the fitted model.

    +
    distimes

    The DT50 and DT90 values for each observed variable.

    +
    SFORB

    If applicable, eigenvalues of SFORB components of the model.

    +

    The print method is called for its side effect, i.e. printing the summary.

    @@ -155,50 +187,366 @@ 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)
-#>  
-#> Error in rxModelVars_(obj): Not compatible with STRSXP: [type=NULL].
-summary(f_saem_dfop_sfo, data = TRUE)
-#> Error in h(simpleError(msg, call)): error in evaluating the argument 'object' in selecting a method for function 'summary': object 'f_saem_dfop_sfo' not found
-# }
-
+    
    # 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)
+print(f_saem_dfop_sfo)
+#> Kinetic nonlinear mixed-effects model fit by SAEM
+#> Structural model:
+#> d_parent/dt = - ((k1 * g * exp(-k1 * time) + k2 * (1 - g) * exp(-k2 *
+#>            time)) / (g * exp(-k1 * time) + (1 - g) * exp(-k2 * time)))
+#>            * parent
+#> d_m1/dt = + f_parent_to_m1 * ((k1 * g * exp(-k1 * time) + k2 * (1 - g)
+#>            * exp(-k2 * time)) / (g * exp(-k1 * time) + (1 - g) *
+#>            exp(-k2 * time))) * parent - k_m1 * m1
+#> 
+#> Data:
+#> 171 observations of 2 variable(s) grouped in 5 datasets
+#> 
+#> Likelihood computed by importance sampling
+#>     AIC   BIC logLik
+#>   828.1 822.7 -400.1
+#> 
+#> Fitted parameters:
+#>                     estimate    lower     upper
+#> parent_0           100.74378 97.81291 103.67465
+#> log_k_m1            -4.06168 -4.17104  -3.95231
+#> f_parent_qlogis     -0.92584 -1.31273  -0.53894
+#> log_k1              -2.81914 -3.60206  -2.03623
+#> log_k2              -3.63916 -4.32672  -2.95161
+#> g_qlogis            -0.02927 -1.15247   1.09394
+#> a.1                  0.86164  0.67928   1.04400
+#> b.1                  0.07973  0.06437   0.09509
+#> SD.parent_0          0.73313 -7.46512   8.93137
+#> SD.log_k_m1          0.06488 -0.06041   0.19017
+#> SD.f_parent_qlogis   0.41955  0.15206   0.68705
+#> SD.log_k1            0.81750  0.29140   1.34361
+#> SD.log_k2            0.75265  0.27939   1.22590
+#> SD.g_qlogis          0.34411 -1.70964   2.39786
+illparms(f_saem_dfop_sfo)
+#> [1] "sd(parent_0)" "sd(log_k_m1)" "sd(g_qlogis)"
+f_saem_dfop_sfo_2 <- update(f_saem_dfop_sfo, covariance.model = diag(c(0, 0, 1, 1, 1, 0)))
+illparms(f_saem_dfop_sfo_2)
+#> character(0)
+intervals(f_saem_dfop_sfo_2)
+#> Approximate 95% confidence intervals
+#> 
+#>  Fixed effects:
+#>                      lower         est.        upper
+#> parent_0       97.54844979 100.46239264 103.37633550
+#> k_m1            0.01575805   0.01729111   0.01897331
+#> f_parent_to_m1  0.21014925   0.28626877   0.37680664
+#> k1              0.02651112   0.05601399   0.11834908
+#> k2              0.01326524   0.02649799   0.05293107
+#> g               0.31467778   0.51297098   0.70726363
+#> 
+#>  Random effects:
+#>                         lower      est.     upper
+#> sd(f_parent_qlogis) 0.1658367 0.4471180 0.7283993
+#> sd(log_k1)          0.2768757 0.7929203 1.3089649
+#> sd(log_k2)          0.2693629 0.7566116 1.2438602
+#> 
+#>  
+#>          lower       est.      upper
+#> a.1 0.70273100 0.88750764 1.07228428
+#> b.1 0.06781347 0.08328016 0.09874685
+summary(f_saem_dfop_sfo_2, data = TRUE)
+#> saemix version used for fitting:      3.1 
+#> mkin version used for pre-fitting:  1.1.2 
+#> R version used for fitting:         4.2.1 
+#> Date of fit:     Wed Aug 10 15:28:15 2022 
+#> Date of summary: Wed Aug 10 15:28:15 2022 
+#> 
+#> Equations:
+#> d_parent/dt = - ((k1 * g * exp(-k1 * time) + k2 * (1 - g) * exp(-k2 *
+#>            time)) / (g * exp(-k1 * time) + (1 - g) * exp(-k2 * time)))
+#>            * parent
+#> d_m1/dt = + f_parent_to_m1 * ((k1 * g * exp(-k1 * time) + k2 * (1 - g)
+#>            * exp(-k2 * time)) / (g * exp(-k1 * time) + (1 - g) *
+#>            exp(-k2 * time))) * parent - k_m1 * m1
+#> 
+#> Data:
+#> 171 observations of 2 variable(s) grouped in 5 datasets
+#> 
+#> Model predictions using solution type analytical 
+#> 
+#> Fitted in 27.674 s
+#> Using 300, 100 iterations and 10 chains
+#> 
+#> Variance model: Two-component variance function 
+#> 
+#> Mean of starting values for individual parameters:
+#>        parent_0        log_k_m1 f_parent_qlogis          log_k1          log_k2 
+#>       101.65645        -4.05368        -0.94311        -2.35943        -4.07006 
+#>        g_qlogis 
+#>        -0.01132 
+#> 
+#> Fixed degradation parameter values:
+#> None
+#> 
+#> Results:
+#> 
+#> Likelihood computed by importance sampling
+#>     AIC   BIC logLik
+#>   825.6 821.3 -401.8
+#> 
+#> Optimised parameters:
+#>                     est.   lower    upper
+#> parent_0        100.4624 97.5484 103.3763
+#> log_k_m1         -4.0576 -4.1504  -3.9647
+#> f_parent_qlogis  -0.9136 -1.3240  -0.5031
+#> log_k1           -2.8822 -3.6302  -2.1341
+#> log_k2           -3.6307 -4.3226  -2.9388
+#> g_qlogis          0.0519 -0.7783   0.8821
+#> 
+#> Correlation: 
+#>                 parnt_0 lg_k_m1 f_prnt_ log_k1  log_k2 
+#> log_k_m1        -0.4102                                
+#> f_parent_qlogis -0.2113  0.2439                        
+#> log_k1           0.1308 -0.1305 -0.0504                
+#> log_k2          -0.0383  0.0592  0.0151  0.0001        
+#> g_qlogis        -0.0029 -0.0118  0.0131 -0.2547 -0.1942
+#> 
+#> Random effects:
+#>                      est.  lower  upper
+#> SD.f_parent_qlogis 0.4471 0.1658 0.7284
+#> SD.log_k1          0.7929 0.2769 1.3090
+#> SD.log_k2          0.7566 0.2694 1.2439
+#> 
+#> Variance model:
+#>        est.   lower   upper
+#> a.1 0.88751 0.70273 1.07228
+#> b.1 0.08328 0.06781 0.09875
+#> 
+#> Backtransformed parameters:
+#>                     est.    lower     upper
+#> parent_0       100.46239 97.54845 103.37634
+#> k_m1             0.01729  0.01576   0.01897
+#> f_parent_to_m1   0.28627  0.21015   0.37681
+#> k1               0.05601  0.02651   0.11835
+#> k2               0.02650  0.01327   0.05293
+#> g                0.51297  0.31468   0.70726
+#> 
+#> Resulting formation fractions:
+#>                 ff
+#> parent_m1   0.2863
+#> parent_sink 0.7137
+#> 
+#> Estimated disappearance times:
+#>         DT50   DT90 DT50back DT50_k1 DT50_k2
+#> parent 17.44  65.15    19.61   12.37   26.16
+#> m1     40.09 133.17       NA      NA      NA
+#> 
+#> Data:
+#>    ds   name time observed  predicted   residual    std standardized
+#>  ds 1 parent    0     89.8  1.005e+02 -10.662393 8.4135    -1.267301
+#>  ds 1 parent    0    104.1  1.005e+02   3.637607 8.4135     0.432355
+#>  ds 1 parent    1     88.7  9.576e+01  -7.063498 8.0244    -0.880249
+#>  ds 1 parent    1     95.5  9.576e+01  -0.263498 8.0244    -0.032837
+#>  ds 1 parent    3     81.8  8.717e+01  -5.369491 7.3135    -0.734185
+#>  ds 1 parent    3     94.5  8.717e+01   7.330509 7.3135     1.002320
+#>  ds 1 parent    7     71.5  7.274e+01  -1.238672 6.1224    -0.202319
+#>  ds 1 parent    7     70.3  7.274e+01  -2.438672 6.1224    -0.398322
+#>  ds 1 parent   14     54.2  5.418e+01   0.022691 4.5984     0.004935
+#>  ds 1 parent   14     49.6  5.418e+01  -4.577309 4.5984    -0.995423
+#>  ds 1 parent   28     31.5  3.241e+01  -0.914545 2.8416    -0.321837
+#>  ds 1 parent   28     28.8  3.241e+01  -3.614545 2.8416    -1.271993
+#>  ds 1 parent   60     12.1  1.283e+01  -0.730904 1.3891    -0.526186
+#>  ds 1 parent   60     13.6  1.283e+01   0.769096 1.3891     0.553681
+#>  ds 1 parent   90      6.2  6.128e+00   0.071981 1.0238     0.070309
+#>  ds 1 parent   90      8.3  6.128e+00   2.171981 1.0238     2.121538
+#>  ds 1 parent  120      2.2  3.022e+00  -0.822164 0.9225    -0.891230
+#>  ds 1 parent  120      2.4  3.022e+00  -0.622164 0.9225    -0.674429
+#>  ds 1     m1    1      0.3  1.163e+00  -0.863423 0.8928    -0.967116
+#>  ds 1     m1    1      0.2  1.163e+00  -0.963423 0.8928    -1.079126
+#>  ds 1     m1    3      2.2  3.233e+00  -1.032930 0.9274    -1.113734
+#>  ds 1     m1    3      3.0  3.233e+00  -0.232930 0.9274    -0.251152
+#>  ds 1     m1    7      6.5  6.495e+00   0.005314 1.0393     0.005113
+#>  ds 1     m1    7      5.0  6.495e+00  -1.494686 1.0393    -1.438116
+#>  ds 1     m1   14     10.2  1.010e+01   0.096372 1.2230     0.078801
+#>  ds 1     m1   14      9.5  1.010e+01  -0.603628 1.2230    -0.493572
+#>  ds 1     m1   28     12.2  1.269e+01  -0.492073 1.3802    -0.356526
+#>  ds 1     m1   28     13.4  1.269e+01   0.707927 1.3802     0.512922
+#>  ds 1     m1   60     11.8  1.086e+01   0.944360 1.2669     0.745420
+#>  ds 1     m1   60     13.2  1.086e+01   2.344360 1.2669     1.850494
+#>  ds 1     m1   90      6.6  7.723e+00  -1.123088 1.0961    -1.024658
+#>  ds 1     m1   90      9.3  7.723e+00   1.576912 1.0961     1.438708
+#>  ds 1     m1  120      3.5  5.184e+00  -1.683936 0.9869    -1.706219
+#>  ds 1     m1  120      5.4  5.184e+00   0.216064 0.9869     0.218923
+#>  ds 2 parent    0    118.0  1.005e+02  17.537607 8.4135     2.084469
+#>  ds 2 parent    0     99.8  1.005e+02  -0.662393 8.4135    -0.078730
+#>  ds 2 parent    1     90.2  9.566e+01  -5.456414 8.0156    -0.680727
+#>  ds 2 parent    1     94.6  9.566e+01  -1.056414 8.0156    -0.131795
+#>  ds 2 parent    3     96.1  8.702e+01   9.082833 7.3009     1.244062
+#>  ds 2 parent    3     78.4  8.702e+01  -8.617167 7.3009    -1.180281
+#>  ds 2 parent    7     77.9  7.298e+01   4.919834 6.1423     0.800981
+#>  ds 2 parent    7     77.7  7.298e+01   4.719834 6.1423     0.768420
+#>  ds 2 parent   14     56.0  5.588e+01   0.124003 4.7372     0.026176
+#>  ds 2 parent   14     54.7  5.588e+01  -1.175997 4.7372    -0.248245
+#>  ds 2 parent   28     36.6  3.719e+01  -0.587869 3.2217    -0.182474
+#>  ds 2 parent   28     36.8  3.719e+01  -0.387869 3.2217    -0.120394
+#>  ds 2 parent   60     22.1  2.013e+01   1.973728 1.8966     1.040673
+#>  ds 2 parent   60     24.7  2.013e+01   4.573728 1.8966     2.411556
+#>  ds 2 parent   90     12.4  1.259e+01  -0.185933 1.3734    -0.135379
+#>  ds 2 parent   90     10.8  1.259e+01  -1.785933 1.3734    -1.300347
+#>  ds 2 parent  120      6.8  7.981e+00  -1.180542 1.1088    -1.064723
+#>  ds 2 parent  120      7.9  7.981e+00  -0.080542 1.1088    -0.072640
+#>  ds 2     m1    1      1.3  1.306e+00  -0.006246 0.8941    -0.006986
+#>  ds 2     m1    3      3.7  3.589e+00   0.110879 0.9365     0.118399
+#>  ds 2     m1    3      4.7  3.589e+00   1.110879 0.9365     1.186217
+#>  ds 2     m1    7      8.1  7.062e+00   1.038045 1.0647     0.974978
+#>  ds 2     m1    7      7.9  7.062e+00   0.838045 1.0647     0.787129
+#>  ds 2     m1   14     10.1  1.065e+01  -0.553713 1.2549    -0.441227
+#>  ds 2     m1   14     10.3  1.065e+01  -0.353713 1.2549    -0.281857
+#>  ds 2     m1   28     10.7  1.284e+01  -2.144854 1.3900    -1.543111
+#>  ds 2     m1   28     12.2  1.284e+01  -0.644854 1.3900    -0.463939
+#>  ds 2     m1   60     10.7  1.082e+01  -0.115278 1.2645    -0.091165
+#>  ds 2     m1   60     12.5  1.082e+01   1.684722 1.2645     1.332337
+#>  ds 2     m1   90      9.1  8.014e+00   1.085607 1.1105     0.977610
+#>  ds 2     m1   90      7.4  8.014e+00  -0.614393 1.1105    -0.553272
+#>  ds 2     m1  120      6.1  5.736e+00   0.363593 1.0079     0.360737
+#>  ds 2     m1  120      4.5  5.736e+00  -1.236407 1.0079    -1.226697
+#>  ds 3 parent    0    106.2  1.005e+02   5.737607 8.4135     0.681955
+#>  ds 3 parent    0    106.9  1.005e+02   6.437607 8.4135     0.765155
+#>  ds 3 parent    1    107.4  9.343e+01  13.972212 7.8311     1.784188
+#>  ds 3 parent    1     96.1  9.343e+01   2.672212 7.8311     0.341229
+#>  ds 3 parent    3     79.4  8.160e+01  -2.196297 6.8531    -0.320484
+#>  ds 3 parent    3     82.6  8.160e+01   1.003703 6.8531     0.146460
+#>  ds 3 parent    7     63.9  6.464e+01  -0.737220 5.4557    -0.135129
+#>  ds 3 parent    7     62.4  6.464e+01  -2.237220 5.4557    -0.410072
+#>  ds 3 parent   14     51.0  4.772e+01   3.278433 4.0722     0.805086
+#>  ds 3 parent   14     47.1  4.772e+01  -0.621567 4.0722    -0.152638
+#>  ds 3 parent   28     36.1  3.303e+01   3.070676 2.8903     1.062400
+#>  ds 3 parent   28     36.6  3.303e+01   3.570676 2.8903     1.235391
+#>  ds 3 parent   60     20.1  1.929e+01   0.808039 1.8355     0.440235
+#>  ds 3 parent   60     19.8  1.929e+01   0.508039 1.8355     0.276789
+#>  ds 3 parent   90     11.3  1.209e+01  -0.794443 1.3425    -0.591785
+#>  ds 3 parent   90     10.7  1.209e+01  -1.394443 1.3425    -1.038728
+#>  ds 3 parent  120      8.2  7.590e+00   0.610002 1.0896     0.559843
+#>  ds 3 parent  120      7.3  7.590e+00  -0.289998 1.0896    -0.266152
+#>  ds 3     m1    0      0.8 -4.263e-14   0.800000 0.8875     0.901401
+#>  ds 3     m1    1      1.8  1.692e+00   0.107665 0.8986     0.119811
+#>  ds 3     m1    1      2.3  1.692e+00   0.607665 0.8986     0.676214
+#>  ds 3     m1    3      4.2  4.455e+00  -0.255347 0.9619    -0.265449
+#>  ds 3     m1    3      4.1  4.455e+00  -0.355347 0.9619    -0.369404
+#>  ds 3     m1    7      6.8  8.124e+00  -1.324338 1.1160    -1.186685
+#>  ds 3     m1    7     10.1  8.124e+00   1.975662 1.1160     1.770309
+#>  ds 3     m1   14     11.4  1.104e+01   0.361860 1.2778     0.283196
+#>  ds 3     m1   14     12.8  1.104e+01   1.761860 1.2778     1.378852
+#>  ds 3     m1   28     11.5  1.177e+01  -0.272554 1.3225    -0.206097
+#>  ds 3     m1   28     10.6  1.177e+01  -1.172554 1.3225    -0.886648
+#>  ds 3     m1   60      7.5  9.242e+00  -1.741667 1.1747    -1.482591
+#>  ds 3     m1   60      8.6  9.242e+00  -0.641667 1.1747    -0.546218
+#>  ds 3     m1   90      7.3  6.837e+00   0.463318 1.0544     0.439398
+#>  ds 3     m1   90      8.1  6.837e+00   1.263318 1.0544     1.198095
+#>  ds 3     m1  120      5.3  4.906e+00   0.394322 0.9770     0.403595
+#>  ds 3     m1  120      3.8  4.906e+00  -1.105678 0.9770    -1.131677
+#>  ds 4 parent    0    104.7  1.005e+02   4.237607 8.4135     0.503670
+#>  ds 4 parent    0     88.3  1.005e+02 -12.162393 8.4135    -1.445587
+#>  ds 4 parent    1     94.2  9.723e+01  -3.029220 8.1458    -0.371877
+#>  ds 4 parent    1     94.6  9.723e+01  -2.629220 8.1458    -0.322772
+#>  ds 4 parent    3     78.1  9.114e+01 -13.041804 7.6420    -1.706592
+#>  ds 4 parent    3     96.5  9.114e+01   5.358196 7.6420     0.701150
+#>  ds 4 parent    7     76.2  8.033e+01  -4.133084 6.7488    -0.612421
+#>  ds 4 parent    7     77.8  8.033e+01  -2.533084 6.7488    -0.375340
+#>  ds 4 parent   14     70.8  6.504e+01   5.757987 5.4889     1.049017
+#>  ds 4 parent   14     67.3  6.504e+01   2.257987 5.4889     0.411371
+#>  ds 4 parent   28     43.1  4.418e+01  -1.080806 3.7849    -0.285557
+#>  ds 4 parent   28     45.1  4.418e+01   0.919194 3.7849     0.242858
+#>  ds 4 parent   60     21.3  2.110e+01   0.200596 1.9686     0.101899
+#>  ds 4 parent   60     23.5  2.110e+01   2.400596 1.9686     1.219459
+#>  ds 4 parent   90     11.8  1.183e+01  -0.034206 1.3263    -0.025791
+#>  ds 4 parent   90     12.1  1.183e+01   0.265794 1.3263     0.200408
+#>  ds 4 parent  120      7.0  6.985e+00   0.014647 1.0612     0.013803
+#>  ds 4 parent  120      6.2  6.985e+00  -0.785353 1.0612    -0.740078
+#>  ds 4     m1    0      1.6 -1.705e-13   1.600000 0.8875     1.802801
+#>  ds 4     m1    1      0.9  6.803e-01   0.219655 0.8893     0.246994
+#>  ds 4     m1    3      3.7  1.927e+00   1.773027 0.9019     1.965880
+#>  ds 4     m1    3      2.0  1.927e+00   0.073027 0.9019     0.080970
+#>  ds 4     m1    7      3.6  4.013e+00  -0.412926 0.9483    -0.435417
+#>  ds 4     m1    7      3.8  4.013e+00  -0.212926 0.9483    -0.224523
+#>  ds 4     m1   14      7.1  6.604e+00   0.495843 1.0441     0.474896
+#>  ds 4     m1   14      6.6  6.604e+00  -0.004157 1.0441    -0.003981
+#>  ds 4     m1   28      9.5  9.077e+00   0.422700 1.1658     0.362576
+#>  ds 4     m1   28      9.3  9.077e+00   0.222700 1.1658     0.191024
+#>  ds 4     m1   60      8.3  8.818e+00  -0.518498 1.1520    -0.450099
+#>  ds 4     m1   60      9.0  8.818e+00   0.181502 1.1520     0.157559
+#>  ds 4     m1   90      6.6  6.738e+00  -0.137785 1.0500    -0.131222
+#>  ds 4     m1   90      7.7  6.738e+00   0.962215 1.0500     0.916383
+#>  ds 4     m1  120      3.7  4.794e+00  -1.093754 0.9732    -1.123914
+#>  ds 4     m1  120      3.5  4.794e+00  -1.293754 0.9732    -1.329429
+#>  ds 5 parent    0    110.4  1.005e+02   9.937607 8.4135     1.181155
+#>  ds 5 parent    0    112.1  1.005e+02  11.637607 8.4135     1.383212
+#>  ds 5 parent    1     93.5  9.372e+01  -0.215694 7.8550    -0.027460
+#>  ds 5 parent    1     91.0  9.372e+01  -2.715694 7.8550    -0.345730
+#>  ds 5 parent    3     71.0  8.226e+01 -11.257156 6.9076    -1.629667
+#>  ds 5 parent    3     89.7  8.226e+01   7.442844 6.9076     1.077480
+#>  ds 5 parent    7     60.4  6.553e+01  -5.128464 5.5289    -0.927571
+#>  ds 5 parent    7     59.1  6.553e+01  -6.428464 5.5289    -1.162699
+#>  ds 5 parent   14     56.5  4.835e+01   8.146351 4.1235     1.975572
+#>  ds 5 parent   14     47.0  4.835e+01  -1.353649 4.1235    -0.328273
+#>  ds 5 parent   28     30.2  3.300e+01  -2.803303 2.8883    -0.970586
+#>  ds 5 parent   28     23.9  3.300e+01  -9.103303 2.8883    -3.151832
+#>  ds 5 parent   60     17.0  1.891e+01  -1.905909 1.8074    -1.054506
+#>  ds 5 parent   60     18.7  1.891e+01  -0.205909 1.8074    -0.113926
+#>  ds 5 parent   90     11.3  1.172e+01  -0.423434 1.3194    -0.320923
+#>  ds 5 parent   90     11.9  1.172e+01   0.176566 1.3194     0.133820
+#>  ds 5 parent  120      9.0  7.281e+00   1.719138 1.0749     1.599402
+#>  ds 5 parent  120      8.1  7.281e+00   0.819138 1.0749     0.762086
+#>  ds 5     m1    0      0.7 -2.842e-13   0.700000 0.8875     0.788726
+#>  ds 5     m1    1      3.0  3.252e+00  -0.252227 0.9279    -0.271821
+#>  ds 5     m1    1      2.6  3.252e+00  -0.652227 0.9279    -0.702895
+#>  ds 5     m1    3      5.1  8.615e+00  -3.515326 1.1413    -3.080237
+#>  ds 5     m1    3      7.5  8.615e+00  -1.115326 1.1413    -0.977283
+#>  ds 5     m1    7     16.5  1.588e+01   0.619041 1.5928     0.388661
+#>  ds 5     m1    7     19.0  1.588e+01   3.119041 1.5928     1.958272
+#>  ds 5     m1   14     22.9  2.189e+01   1.014705 2.0272     0.500543
+#>  ds 5     m1   14     23.2  2.189e+01   1.314705 2.0272     0.648529
+#>  ds 5     m1   28     22.2  2.369e+01  -1.487604 2.1632    -0.687701
+#>  ds 5     m1   28     24.4  2.369e+01   0.712396 2.1632     0.329332
+#>  ds 5     m1   60     15.5  1.869e+01  -3.193942 1.7920    -1.782295
+#>  ds 5     m1   60     19.8  1.869e+01   1.106058 1.7920     0.617206
+#>  ds 5     m1   90     14.9  1.380e+01   1.103454 1.4518     0.760041
+#>  ds 5     m1   90     14.2  1.380e+01   0.403454 1.4518     0.277892
+#>  ds 5     m1  120     10.9  9.864e+00   1.035963 1.2093     0.856637
+#>  ds 5     m1  120     10.4  9.864e+00   0.535963 1.2093     0.443187
+# }
+
    @@ -213,7 +561,7 @@ saemix authors for the parts inherited from saemix.

    -

    Site built with pkgdown 2.0.2.

    +

    Site built with pkgdown 2.0.6.

    -- cgit v1.2.1