Estimation of parameter distributions from mmkin row objects

Date: Wed, 15 Apr 2020 01:27:02 +0200 Subject: A plot method for nlme.mmkin fits --- docs/reference/index.html | 23 ++- docs/reference/nlme.html | 32 ++-- docs/reference/nlme.mmkin.html | 301 +++++++++++++++++++++++++++++++++++ docs/reference/plot.nlme.mmkin-1.png | Bin 0 -> 32297 bytes docs/reference/plot.nlme.mmkin.html | 252 +++++++++++++++++++++++++++++ 5 files changed, 592 insertions(+), 16 deletions(-) create mode 100644 docs/reference/nlme.mmkin.html create mode 100644 docs/reference/plot.nlme.mmkin-1.png create mode 100644 docs/reference/plot.nlme.mmkin.html (limited to 'docs/reference') diff --git a/docs/reference/index.html b/docs/reference/index.html index 31f7678c..3c5f1b38 100644 --- a/docs/reference/index.html +++ b/docs/reference/index.html @@ -270,12 +270,32 @@ of an mmkin object

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

Calculate the AIC for a column of an mmkin object

+ + + + +

Mixed models

Create and work with nonlinear mixed models

+ + + + + +

nlme(<mmkin>)

+ +

Create an nlme model for an mmkin row object

+ + + +

plot(<nlme.mmkin>)

+ +

Plot a fitted nonlinear mixed model obtained via an mmkin row object

nlme_function() mean_degparms() nlme_data()

Estimation of parameter distributions from mmkin row objects

Helper functions to create nlme models from mmkin row objects

@@ -576,6 +596,7 @@ kinetic models fitted with mkinfit

Main functions

Show results

Work with mmkin objects

Mixed models

Datasets and known results

NAFTA guidance

Helper functions mainly used internally

diff --git a/docs/reference/nlme.html b/docs/reference/nlme.html index 1b05f882..981845fe 100644 --- a/docs/reference/nlme.html +++ b/docs/reference/nlme.html @@ -6,7 +6,7 @@ -Estimation of parameter distributions from mmkin row objects — nlme_function • mkin +Helper functions to create nlme models from mmkin row objects — nlme_function • mkin @@ -35,7 +35,7 @@ - +

@@ -255,37 +255,39 @@ datasets.

nlme_f_sfo_sfo <- nlme_function(f_2["SFO-SFO", ]) nlme_f_sfo_sfo_ff <- nlme_function(f_2["SFO-SFO-ff", ]) nlme_f_fomc_sfo <- nlme_function(f_2["FOMC-SFO", ]) + assign("nlme_f_sfo_sfo", nlme_f_sfo_sfo, globalenv()) + assign("nlme_f_sfo_sfo_ff", nlme_f_sfo_sfo_ff, globalenv()) + assign("nlme_f_fomc_sfo", nlme_f_fomc_sfo, globalenv()) - # Allowing for correlations between random effects leads to non-convergence + # Allowing for correlations between random effects (not shown) + # leads to non-convergence f_nlme_sfo_sfo <- nlme(value ~ nlme_f_sfo_sfo(name, time, parent_0, log_k_parent_sink, log_k_parent_A1, log_k_A1_sink), data = grouped_data_2, fixed = parent_0 + log_k_parent_sink + log_k_parent_A1 + log_k_A1_sink ~ 1, random = pdDiag(parent_0 + log_k_parent_sink + log_k_parent_A1 + log_k_A1_sink ~ 1), - start = mean_dp_sfo_sfo)

#> Error in nlme_f_sfo_sfo(name, time, parent_0, log_k_parent_sink, log_k_parent_A1, log_k_A1_sink): konnte Funktion "nlme_f_sfo_sfo" nicht finden

+ start = mean_dp_sfo_sfo) + # augPred does not see to work on this object, so no plot is shown + # The same model fitted with transformed formation fractions does not converge f_nlme_sfo_sfo_ff <- nlme(value ~ nlme_f_sfo_sfo_ff(name, time, parent_0, log_k_parent, log_k_A1, f_parent_ilr_1), data = grouped_data_2, fixed = parent_0 + log_k_parent + log_k_A1 + f_parent_ilr_1 ~ 1, random = pdDiag(parent_0 + log_k_parent + log_k_A1 + f_parent_ilr_1 ~ 1), - start = mean_dp_sfo_sfo_ff)

#> Error in nlme_f_sfo_sfo_ff(name, time, parent_0, log_k_parent, log_k_A1, f_parent_ilr_1): konnte Funktion "nlme_f_sfo_sfo_ff" nicht finden

- # It does converge with this version of reduced random effects - f_nlme_sfo_sfo_ff <- nlme(value ~ nlme_f_sfo_sfo_ff(name, time, - parent_0, log_k_parent, log_k_A1, f_parent_ilr_1), - data = grouped_data_2, - fixed = parent_0 + log_k_parent + log_k_A1 + f_parent_ilr_1 ~ 1, - random = pdDiag(parent_0 + log_k_parent ~ 1), - start = mean_dp_sfo_sfo_ff)

#> Error in nlme_f_sfo_sfo_ff(name, time, parent_0, log_k_parent, log_k_A1, f_parent_ilr_1): konnte Funktion "nlme_f_sfo_sfo_ff" nicht finden

+ start = mean_dp_sfo_sfo_ff)

#> Error in nlme.formula(value ~ nlme_f_sfo_sfo_ff(name, time, parent_0, log_k_parent, log_k_A1, f_parent_ilr_1), data = grouped_data_2, fixed = parent_0 + log_k_parent + log_k_A1 + f_parent_ilr_1 ~ 1, random = pdDiag(parent_0 + log_k_parent + log_k_A1 + f_parent_ilr_1 ~ 1), start = mean_dp_sfo_sfo_ff): step halving factor reduced below minimum in PNLS step

f_nlme_fomc_sfo <- nlme(value ~ nlme_f_fomc_sfo(name, time, parent_0, log_alpha, log_beta, log_k_A1, f_parent_ilr_1), data = grouped_data_2, fixed = parent_0 + log_alpha + log_beta + log_k_A1 + f_parent_ilr_1 ~ 1, random = pdDiag(parent_0 + log_alpha + log_beta + log_k_A1 + f_parent_ilr_1 ~ 1), - start = mean_dp_fomc_sfo)

#> Error in nlme_f_fomc_sfo(name, time, parent_0, log_alpha, log_beta, log_k_A1, f_parent_ilr_1): konnte Funktion "nlme_f_fomc_sfo" nicht finden

+ start = mean_dp_fomc_sfo) + # DFOP-SFO and SFORB-SFO did not converge with full random effects - anova(f_nlme_fomc_sfo, f_nlme_sfo_sfo)

#> Error in anova(f_nlme_fomc_sfo, f_nlme_sfo_sfo): Objekt 'f_nlme_fomc_sfo' nicht gefunden

# } + anova(f_nlme_fomc_sfo, f_nlme_sfo_sfo)

#> Model df AIC BIC logLik Test L.Ratio p-value +#> f_nlme_fomc_sfo 1 11 932.5817 967.0755 -455.2909 +#> f_nlme_sfo_sfo 2 9 1089.2492 1117.4714 -535.6246 1 vs 2 160.6675 <.0001

# }

diff --git a/docs/reference/nlme.mmkin.html b/docs/reference/nlme.mmkin.html new file mode 100644 index 00000000..e1b1ff77 --- /dev/null +++ b/docs/reference/nlme.mmkin.html @@ -0,0 +1,301 @@ + + + + + + + + +Create an nlme model for an mmkin row object — nlme.mmkin • mkin + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +

+ + + + +

+ +

+ + +

This functions sets up a nonlinear mixed effects model for an mmkin row +object. An mmkin row object is essentially a list of mkinfit objects that +have been obtained by fitting the same model to a list of datasets.

+ +

# S3 method for mmkin
+nlme(
+  model,
+  data = sys.frame(sys.parent()),
+  fixed,
+  random = fixed,
+  groups,
+  start,
+  correlation = NULL,
+  weights = NULL,
+  subset,
+  method = c("ML", "REML"),
+  na.action = na.fail,
+  naPattern,
+  control = list(),
+  verbose = FALSE
+)

+ +

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, all fixed effects are complemented +with uncorrelated random effects
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

+ +

Value

+ +

Upon success, a fitted nlme.mmkin object, which is an nlme object + with additional elements

Examples

ds <- lapply(experimental_data_for_UBA_2019[6:10],
+ function(x) subset(x$data[c("name", "time", "value")], name == "parent"))
+f <- mmkin("SFO", ds, quiet = TRUE, cores = 1)
+library(nlme)
+f_nlme <- nlme(f)
+nlme(f, random = parent_0 ~ 1)
#> Nonlinear mixed-effects model fit by maximum likelihood
+#>   Model: value ~ deg_func(name, time, parent_0, log_k_parent_sink) 
+#>   Data: structure(list(ds = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,  1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L,  2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L,  3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L,  4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L,  4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L,  5L, 5L), .Label = c("1", "2", "3", "4", "5"), class = c("ordered",  "factor")), name = c("parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent"), time = c(0, 0, 3, 3, 6, 6, 10, 10, 20, 20,  34, 34, 55, 55, 90, 90, 112, 112, 132, 132, 0, 0, 3, 3, 7, 7,  14, 14, 30, 30, 60, 60, 90, 90, 120, 120, 180, 180, 0, 0, 1,  1, 3, 3, 8, 8, 14, 14, 27, 27, 48, 48, 70, 70, 0, 0, 1, 1, 3,  3, 8, 8, 14, 14, 27, 27, 48, 48, 70, 70, 91, 91, 120, 120, 0,  0, 8, 8, 14, 14, 21, 21, 41, 41, 63, 63, 91, 91, 120, 120), value = c(97.2,  96.4, 71.1, 69.2, 58.1, 56.6, 44.4, 43.4, 33.3, 29.2, 17.6, 18,  10.5, 9.3, 4.5, 4.7, 3, 3.4, 2.3, 2.7, 93.6, 92.3, 87, 82.2,  74, 73.9, 64.2, 69.5, 54, 54.6, 41.1, 38.4, 32.5, 35.5, 28.1,  29, 26.5, 27.6, 91.9, 90.8, 64.9, 66.2, 43.5, 44.1, 18.3, 18.1,  10.2, 10.8, 4.9, 3.3, 1.6, 1.5, 1.1, 0.9, 99.8, 98.3, 77.1, 77.2,  59, 58.1, 27.4, 29.2, 19.1, 29.6, 10.1, 18.2, 4.5, 9.1, 2.3,  2.9, 2, 1.8, 2, 2.2, 96.1, 94.3, 73.9, 73.9, 69.4, 73.1, 65.6,  65.3, 55.9, 54.4, 47, 49.3, 44.7, 46.7, 42.1, 41.3)), row.names = c(NA,  -90L), class = c("nfnGroupedData", "nfGroupedData", "groupedData",  "data.frame"), formula = value ~ time | ds, FUN = function (x)  max(x, na.rm = TRUE), order.groups = FALSE) 
+#>   Log-likelihood: -394.4901
+#>   Fixed: list(parent_0 ~ 1, log_k_parent_sink ~ 1) 
+#>          parent_0 log_k_parent_sink 
+#>         73.985522         -3.869079 
+#> 
+#> Random effects:
+#>  Formula: parent_0 ~ 1 | ds
+#>         parent_0 Residual
+#> StdDev:  18.6134 18.22029
+#> 
+#> Number of Observations: 90
+#> Number of Groups: 5 
f_nlme <- nlme(f, start = c(parent_0 = 100, log_k_parent_sink = 0.1))
+update(f_nlme, random = parent_0 ~ 1)
#> Nonlinear mixed-effects model fit by maximum likelihood
+#>   Model: value ~ deg_func(name, time, parent_0, log_k_parent_sink) 
+#>   Data: structure(list(ds = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,  1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L,  2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L,  3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L,  4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L,  4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L,  5L, 5L), .Label = c("1", "2", "3", "4", "5"), class = c("ordered",  "factor")), name = c("parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent", "parent", "parent", "parent", "parent", "parent",  "parent", "parent"), time = c(0, 0, 3, 3, 6, 6, 10, 10, 20, 20,  34, 34, 55, 55, 90, 90, 112, 112, 132, 132, 0, 0, 3, 3, 7, 7,  14, 14, 30, 30, 60, 60, 90, 90, 120, 120, 180, 180, 0, 0, 1,  1, 3, 3, 8, 8, 14, 14, 27, 27, 48, 48, 70, 70, 0, 0, 1, 1, 3,  3, 8, 8, 14, 14, 27, 27, 48, 48, 70, 70, 91, 91, 120, 120, 0,  0, 8, 8, 14, 14, 21, 21, 41, 41, 63, 63, 91, 91, 120, 120), value = c(97.2,  96.4, 71.1, 69.2, 58.1, 56.6, 44.4, 43.4, 33.3, 29.2, 17.6, 18,  10.5, 9.3, 4.5, 4.7, 3, 3.4, 2.3, 2.7, 93.6, 92.3, 87, 82.2,  74, 73.9, 64.2, 69.5, 54, 54.6, 41.1, 38.4, 32.5, 35.5, 28.1,  29, 26.5, 27.6, 91.9, 90.8, 64.9, 66.2, 43.5, 44.1, 18.3, 18.1,  10.2, 10.8, 4.9, 3.3, 1.6, 1.5, 1.1, 0.9, 99.8, 98.3, 77.1, 77.2,  59, 58.1, 27.4, 29.2, 19.1, 29.6, 10.1, 18.2, 4.5, 9.1, 2.3,  2.9, 2, 1.8, 2, 2.2, 96.1, 94.3, 73.9, 73.9, 69.4, 73.1, 65.6,  65.3, 55.9, 54.4, 47, 49.3, 44.7, 46.7, 42.1, 41.3)), row.names = c(NA,  -90L), class = c("nfnGroupedData", "nfGroupedData", "groupedData",  "data.frame"), formula = value ~ time | ds, FUN = function (x)  max(x, na.rm = TRUE), order.groups = FALSE) 
+#>   Log-likelihood: -404.3729
+#>   Fixed: list(parent_0 ~ 1, log_k_parent_sink ~ 1) 
+#>          parent_0 log_k_parent_sink 
+#>         75.933480         -3.555983 
+#> 
+#> Random effects:
+#>  Formula: parent_0 ~ 1 | ds
+#>            parent_0 Residual
+#> StdDev: 0.002416792 21.63027
+#> 
+#> Number of Observations: 90
+#> Number of Groups: 5

+ +

+ + + +

+ + + + + + + + diff --git a/docs/reference/plot.nlme.mmkin-1.png b/docs/reference/plot.nlme.mmkin-1.png new file mode 100644 index 00000000..0717f30d Binary files /dev/null and b/docs/reference/plot.nlme.mmkin-1.png differ diff --git a/docs/reference/plot.nlme.mmkin.html b/docs/reference/plot.nlme.mmkin.html new file mode 100644 index 00000000..d5b7c00c --- /dev/null +++ b/docs/reference/plot.nlme.mmkin.html @@ -0,0 +1,252 @@ + + + + + + + + +Plot a fitted nonlinear mixed model obtained via an mmkin row object — plot.nlme.mmkin • mkin + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +

+ + + + +

+ +

+ + +

Plot a fitted nonlinear mixed model obtained via an mmkin row object

+ +

# S3 method for nlme.mmkin
+plot(
+  x,
+  i = 1:ncol(x$mmkin_orig),
+  main = "auto",
+  legends = 1,
+  resplot = c("time", "errmod"),
+  standardized = FALSE,
+  cex = 0.7,
+  rel.height.middle = 0.9,
+  ymax = "auto",
+  ...
+)

+ +

Arguments

+ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +

x	An object of class `nlme.mmkin`
i	A numeric index to select datasets for which to plot the nlme fit, +in case plots get too large
main	The main title placed on the outer margin of the plot.
legends	An index for the fits for which legends should be shown.
resplot	Should the residuals plotted against time, using +`mkinresplot`, or as squared residuals against predicted +values, with the error model, using `mkinerrplot`.
standardized	Should the residuals be standardized? This option +is passed to `mkinresplot`, it only takes effect if +`resplot = "time"`.
cex	Passed to the plot functions and `mtext`.
rel.height.middle	The relative height of the middle plot, if more +than two rows of plots are shown.
ymax	Maximum y axis value for `plot.mkinfit`.
...	Further arguments passed to `plot.mkinfit` and +`mkinresplot`.

+ +

Value

+ +

The function is called for its side effect.

+ +

Examples

ds <- lapply(experimental_data_for_UBA_2019[6:10],
+ function(x) subset(x$data[c("name", "time", "value")], name == "parent"))
+f <- mmkin("SFO", ds, quiet = TRUE, cores = 1)
+#plot(f) # too many panels for pkgdown
+library(nlme)
+f_nlme <- nlme(f)
+
+#plot(f_nlme) # too many panels for pkgdown
+plot(f_nlme, 1:2)

+ +

+ + + +

+ + + + + + + + -- cgit v1.2.1

Mixed models

Estimation of parameter distributions from mmkin row objects

Helper functions to create nlme models from mmkin row objects

Create an nlme model for an mmkin row object

Arguments

Value

See also

Examples

Contents

Plot a fitted nonlinear mixed model obtained via an mmkin row object

Arguments

Value

Examples

Contents

Author