Create degradation functions for known analytical solutions

From 92bd33824bde6b6b21bfc7e30953092a74d3cce5 Mon Sep 17 00:00:00 2001 From: Johannes Ranke Date: Thu, 7 May 2020 22:13:33 +0200 Subject: Another overhaul of analytical solutions Still in preparation for analytical solutions of coupled models --- docs/reference/DFOP.solution.html | 2 +- docs/reference/FOMC.solution.html | 2 +- docs/reference/HS.solution.html | 33 +++++- docs/reference/SFO.solution.html | 2 +- docs/reference/SFORB.solution.html | 22 +++- docs/reference/create_deg_func.html | 198 ++++++++++++++++++++++++++++++++++ docs/reference/index.html | 4 +- docs/reference/logistic.solution.html | 63 ++--------- docs/reference/mkinmod.html | 16 +-- docs/reference/mkinpredict.html | 48 ++++++--- docs/reference/nlme-1.png | Bin 70555 -> 71651 bytes docs/reference/nlme.html | 75 ++++++------- 12 files changed, 347 insertions(+), 118 deletions(-) create mode 100644 docs/reference/create_deg_func.html (limited to 'docs/reference') diff --git a/docs/reference/DFOP.solution.html b/docs/reference/DFOP.solution.html index 46e0cc91..84e33f9e 100644 --- a/docs/reference/DFOP.solution.html +++ b/docs/reference/DFOP.solution.html @@ -73,7 +73,7 @@ two exponential decline functions." /> diff --git a/docs/reference/FOMC.solution.html b/docs/reference/FOMC.solution.html index 6c1f5229..369398d5 100644 --- a/docs/reference/FOMC.solution.html +++ b/docs/reference/FOMC.solution.html @@ -73,7 +73,7 @@ a decreasing rate constant." /> diff --git a/docs/reference/HS.solution.html b/docs/reference/HS.solution.html index ad33d948..2684ee60 100644 --- a/docs/reference/HS.solution.html +++ b/docs/reference/HS.solution.html @@ -73,7 +73,7 @@ between them." /> @@ -151,6 +151,22 @@ between them.

Arguments

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

t	Time.
parent_0	Starting value for the response variable at time zero.
k1	First kinetic constant.
k2	Second kinetic constant.
tb	Break point. Before this time, exponential decline according to @@ -159,6 +175,21 @@ according to `k2`.

Value

+ +

The value of the response variable at time t.

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, + EC Document Reference Sanco/10058/2005 version 2.0, 434 pp, + http://esdac.jrc.ec.europa.eu/projects/degradation-kinetics +FOCUS (2014) “Generic guidance for Estimating Persistence + and Degradation Kinetics from Environmental Fate Studies on Pesticides in + EU Registration” Report of the FOCUS Work Group on Degradation Kinetics, + Version 1.1, 18 December 2014 + http://esdac.jrc.ec.europa.eu/projects/degradation-kinetics

Arguments

+ + + + + + + + @@ -176,6 +184,18 @@ fraction.

The value of the response variable, which is the sum of free and bound fractions at time t.

References

+ +

Arguments

t	Time.
parent_0	Starting value for the response variable at time zero.
k_12	Kinetic constant describing transfer from free to bound.

+ + + + + + + + + +

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"))
#> Successfully compiled differential equation model from auto-generated C code.
fit <- mkinfit(SFO_SFO, FOCUS_2006_D, quiet = TRUE)
#> Warning: Observations with value of zero were removed from the data

+ + + + + + + + + + + + + + + diff --git a/docs/reference/index.html b/docs/reference/index.html index f77881fe..75820166 100644 --- a/docs/reference/index.html +++ b/docs/reference/index.html @@ -71,7 +71,7 @@ @@ -609,7 +609,7 @@ kinetic models fitted with mkinfit

logistic.solution() logistic.solution()

logistic.solution()

Logistic kinetics

diff --git a/docs/reference/logistic.solution.html b/docs/reference/logistic.solution.html index 9697d420..b0ec77b1 100644 --- a/docs/reference/logistic.solution.html +++ b/docs/reference/logistic.solution.html @@ -41,8 +41,6 @@ @@ -75,7 +73,7 @@ an increasing rate constant, supposedly caused by microbial growth" /> @@ -139,20 +137,16 @@ an increasing rate constant, supposedly caused by microbial growth" />

Function describing exponential decline from a defined starting value, with -an increasing rate constant, supposedly caused by microbial growth

Function describing exponential decline from a defined starting value, with an increasing rate constant, supposedly caused by microbial growth

logistic.solution(t, parent_0, kmax, k0, r)
-
-logistic.solution(t, parent_0, kmax, k0, r)

logistic.solution(t, parent_0, kmax, k0, r)

Arguments

@@ -177,10 +171,6 @@ an increasing rate constant, supposedly caused by microbial growth

- - - -

r	Growth rate of the increase in the rate constant.
parent.0	Starting value for the response variable at time zero.

Value

@@ -189,12 +179,15 @@ an increasing rate constant, supposedly caused by microbial growth

Note

The solution of the logistic model reduces to the - SFO.solution if k0 is equal to kmax.

The solution of the logistic model reduces to the SFO.solution if k0 is equal to kmax.

References

FOCUS (2014) “Generic guidance for Estimating Persistence +

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, + EC Document Reference Sanco/10058/2005 version 2.0, 434 pp, + http://esdac.jrc.ec.europa.eu/projects/degradation-kinetics +FOCUS (2014) “Generic guidance for Estimating Persistence and Degradation Kinetics from Environmental Fate Studies on Pesticides in EU Registration” Report of the FOCUS Work Group on Degradation Kinetics, Version 1.1, 18 December 2014 @@ -247,44 +240,6 @@ an increasing rate constant, supposedly caused by microbial growth

#> k0 0.4448749 #> r 1.1821120 #> sigma 7.3256566

endpoints(m)$distimes

#> DT50 DT90 DT50_k0 DT50_kmax -#> parent 36.86533 62.41511 4297.853 10.83349

- - # Reproduce the plot on page 57 of FOCUS (2014) - plot(function(x) logistic.solution(x, 100, 0.08, 0.0001, 0.2), - from = 0, to = 100, ylim = c(0, 100), - xlab = "Time", ylab = "Residue")

plot(function(x) logistic.solution(x, 100, 0.08, 0.0001, 0.4), - from = 0, to = 100, add = TRUE, lty = 2, col = 2)

plot(function(x) logistic.solution(x, 100, 0.08, 0.0001, 0.8), - from = 0, to = 100, add = TRUE, lty = 3, col = 3)

plot(function(x) logistic.solution(x, 100, 0.08, 0.001, 0.2), - from = 0, to = 100, add = TRUE, lty = 4, col = 4)

plot(function(x) logistic.solution(x, 100, 0.08, 0.08, 0.2), - from = 0, to = 100, add = TRUE, lty = 5, col = 5)

legend("topright", inset = 0.05, - legend = paste0("k0 = ", c(0.0001, 0.0001, 0.0001, 0.001, 0.08), - ", r = ", c(0.2, 0.4, 0.8, 0.2, 0.2)), - lty = 1:5, col = 1:5)

- # Fit with synthetic data - logistic <- mkinmod(parent = mkinsub("logistic")) - - sampling_times = c(0, 1, 3, 7, 14, 28, 60, 90, 120) - parms_logistic <- c(kmax = 0.08, k0 = 0.0001, r = 0.2) - d_logistic <- mkinpredict(logistic, - parms_logistic, c(parent = 100), - sampling_times) - d_2_1 <- add_err(d_logistic, - sdfunc = function(x) sigma_twocomp(x, 0.5, 0.07), - n = 1, reps = 2, digits = 5, LOD = 0.1, seed = 123456)[[1]] - - m <- mkinfit("logistic", d_2_1, quiet = TRUE) - plot_sep(m)

summary(m)$bpar

#> Estimate se_notrans t value Pr(>t) Lower -#> parent_0 1.057896e+02 1.9023449703 55.610119 3.768361e-16 1.016451e+02 -#> kmax 6.398190e-02 0.0143201031 4.467978 3.841829e-04 3.929235e-02 -#> k0 1.612775e-04 0.0005866813 0.274898 3.940351e-01 5.846688e-08 -#> r 2.263946e-01 0.1718110715 1.317695 1.061044e-01 4.335843e-02 -#> sigma 5.332935e+00 0.9145907310 5.830952 4.036926e-05 3.340213e+00 -#> Upper -#> parent_0 109.9341588 -#> kmax 0.1041853 -#> k0 0.4448749 -#> r 1.1821120 -#> sigma 7.3256566

endpoints(m)$distimes

#> DT50 DT90 DT50_k0 DT50_kmax #> parent 36.86533 62.41511 4297.853 10.83349

diff --git a/docs/reference/mkinmod.html b/docs/reference/mkinmod.html index d1204e2e..2f2e89d9 100644 --- a/docs/reference/mkinmod.html +++ b/docs/reference/mkinmod.html @@ -75,7 +75,7 @@ list of lists can be given in the speclist argument." /> @@ -152,7 +152,7 @@ list of lists can be given in the speclist argument.

mkinmod(
   ...,
-  use_of_ff = "min",
+  use_of_ff = "max",
   speclist = NULL,
   quiet = FALSE,
   verbose = FALSE
@@ -252,15 +252,15 @@ in the FOCUS and NAFTA guidance documents are used.
 SFO_SFO <- mkinmod(
   parent = mkinsub("SFO", "m1"),
   m1 = mkinsub("SFO"), verbose = TRUE)#> Compilation argument:
-#>  /usr/lib/R/bin/R CMD SHLIB file66a9718e919b.c 2> file66a9718e919b.c.err.txt 
+#>  /usr/lib/R/bin/R CMD SHLIB fileb6a4eaab60.c 2> fileb6a4eaab60.c.err.txt 
 #> Program source:
 #>   1: #include <R.h>
 #>   2: 
 #>   3: 
 #>   4: static double parms [3];
-#>   5: #define k_parent_sink parms[0]
-#>   6: #define k_parent_m1 parms[1]
-#>   7: #define k_m1_sink parms[2]
+#>   5: #define k_parent parms[0]
+#>   6: #define f_parent_to_m1 parms[1]
+#>   7: #define k_m1 parms[2]
 #>   8: 
 #>   9: void initpar(void (* odeparms)(int *, double *)) {
 #>  10:     int N = 3;
@@ -270,8 +270,8 @@ in the FOCUS and NAFTA guidance documents are used.
 #>  14: 
 #>  15: void func ( int * n, double * t, double * y, double * f, double * rpar, int * ipar ) {
 #>  16: 
-#>  17: f[0] = - k_parent_sink * y[0] - k_parent_m1 * y[0];
-#>  18: f[1] = + k_parent_m1 * y[0] - k_m1_sink * y[1];
+#>  17: f[0] = - k_parent * y[0];
+#>  18: f[1] = + f_parent_to_m1 * k_parent * y[0] - k_m1 * y[1];
 #>  19: }
#> Successfully compiled differential equation model from auto-generated C code.

 # If we have several parallel metabolites
 # (compare tests/testthat/test_synthetic_data_for_UBA_2014.R)
diff --git a/docs/reference/mkinpredict.html b/docs/reference/mkinpredict.html
index 689fb7c7..21c13156 100644
--- a/docs/reference/mkinpredict.html
+++ b/docs/reference/mkinpredict.html
@@ -74,7 +74,7 @@ kinetic parameters and initial values for the state variables." />
       
       
         mkin
-        0.9.49.11
+        0.9.50
       
     
 
@@ -268,8 +268,29 @@ solver is used.
     
 SFO <- mkinmod(degradinol = mkinsub("SFO"))
 # Compare solution types
-mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20,
-      solution_type = "analytical")
#> Error in (function (t, parent_0, k) {    parent = parent_0 * exp(-k * t)})(t = 0:20, parent.0 = c(degradinol = 100), k = 0.3): unbenutztes Argument (parent.0 = 100)
mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20,
+mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100), 0:20,
+      solution_type = "analytical")
#>    time  degradinol
+#> 1     0 100.0000000
+#> 2     1  74.0818221
+#> 3     2  54.8811636
+#> 4     3  40.6569660
+#> 5     4  30.1194212
+#> 6     5  22.3130160
+#> 7     6  16.5298888
+#> 8     7  12.2456428
+#> 9     8   9.0717953
+#> 10    9   6.7205513
+#> 11   10   4.9787068
+#> 12   11   3.6883167
+#> 13   12   2.7323722
+#> 14   13   2.0241911
+#> 15   14   1.4995577
+#> 16   15   1.1108997
+#> 17   16   0.8229747
+#> 18   17   0.6096747
+#> 19   18   0.4516581
+#> 20   19   0.3345965
+#> 21   20   0.2478752
mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100), 0:20,
       solution_type = "deSolve")
#>    time  degradinol
 #> 1     0 100.0000000
 #> 2     1  74.0818221
@@ -291,7 +312,7 @@ solver is used.
 #> 18   17   0.6096747
 #> 19   18   0.4516581
 #> 20   19   0.3345965
-#> 21   20   0.2478752
mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20,
+#> 21   20   0.2478752
mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100), 0:20,
       solution_type = "deSolve", use_compiled = FALSE)
#>    time  degradinol
 #> 1     0 100.0000000
 #> 2     1  74.0818221
@@ -313,7 +334,7 @@ solver is used.
 #> 18   17   0.6096747
 #> 19   18   0.4516581
 #> 20   19   0.3345965
-#> 21   20   0.2478752
mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20,
+#> 21   20   0.2478752
mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100), 0:20,
       solution_type = "eigen")
#>    time  degradinol
 #> 1     0 100.0000000
 #> 2     1  74.0818221
@@ -337,25 +358,26 @@ solver is used.
 #> 20   19   0.3345965
 #> 21   20   0.2478752

 # Compare integration methods to analytical solution
-mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20,
-      solution_type = "analytical")[21,]
#> Error in (function (t, parent_0, k) {    parent = parent_0 * exp(-k * t)})(t = 0:20, parent.0 = c(degradinol = 100), k = 0.3): unbenutztes Argument (parent.0 = 100)
mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20,
+mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100), 0:20,
+      solution_type = "analytical")[21,]
#>    time degradinol
+#> 21   20  0.2478752
mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100), 0:20,
       method = "lsoda")[21,]
#>    time degradinol
-#> 21   20  0.2478752
mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20,
+#> 21   20  0.2478752
mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100), 0:20,
       method = "ode45")[21,]
#>    time degradinol
-#> 21   20  0.2478752
mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20,
+#> 21   20  0.2478752
mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100), 0:20,
       method = "rk4")[21,]
#>    time degradinol
 #> 21   20  0.2480043
# rk4 is not as precise here
 
 # The number of output times used to make a lot of difference until the
 # default for atol was adjusted
-mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100),
+mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100),
       seq(0, 20, by = 0.1))[201,]
#>     time degradinol
-#> 201   20  0.2478752
mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100),
+#> 201   20  0.2478752
mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100),
       seq(0, 20, by = 0.01))[2001,]
#>      time degradinol
 #> 2001   20  0.2478752

 # Check compiled model versions - they are faster than the eigenvalue based solutions!
 SFO_SFO = mkinmod(parent = list(type = "SFO", to = "m1"),
-                  m1 = list(type = "SFO"))
#> Successfully compiled differential equation model from auto-generated C code.
if(require(rbenchmark)) {
+                  m1 = list(type = "SFO"), use_of_ff = "min")
#> Successfully compiled differential equation model from auto-generated C code.
if(require(rbenchmark)) {
   benchmark(
     eigen = mkinpredict(SFO_SFO, c(k_parent_m1 = 0.05, k_parent_sink = 0.1, k_m1_sink = 0.01),
       c(parent = 100, m1 = 0), seq(0, 20, by = 0.1),
@@ -371,7 +393,7 @@ solver is used.
 }
#> Lade nötiges Paket: rbenchmark
#>               test replications elapsed relative user.self sys.self user.child
 #> 3          deSolve           10   0.229   28.625     0.229        0          0
 #> 2 deSolve_compiled           10   0.008    1.000     0.008        0          0
-#> 1            eigen           10   0.025    3.125     0.026        0          0
+#> 1            eigen           10   0.026    3.250     0.025        0          0
 #>   sys.child
 #> 3         0
 #> 2         0
diff --git a/docs/reference/nlme-1.png b/docs/reference/nlme-1.png
index 68ccb43f..8db1f999 100644
Binary files a/docs/reference/nlme-1.png and b/docs/reference/nlme-1.png differ
diff --git a/docs/reference/nlme.html b/docs/reference/nlme.html
index 70c6b63c..b92d2141 100644
--- a/docs/reference/nlme.html
+++ b/docs/reference/nlme.html
@@ -10,23 +10,27 @@
 
 
 
-
+
 
 
-
+
 
-
+
+
+
+
+
 
 
-
-
+
+
 
 
-
+
 
 
-
-
+
+
 
 
 
@@ -40,7 +44,6 @@
 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." />
-
 
 
 
@@ -58,7 +61,7 @@ datasets." />
 
   
 
-  
+  
     
       
       
@@ -72,7 +75,7 @@ datasets." />
       
       
         mkin
-        0.9.49.10
+        0.9.50
       
     
 
@@ -116,7 +119,12 @@ datasets." />
 
       
       
-        
+        
+  
+    
+     
+  
+
       
       
     

@@ -131,7 +139,7 @@ datasets." />
   
     
     Helper functions to create nlme models from mmkin row objects
-    
+    Source: R/nlme.R
     nlme.Rd
     
 
@@ -177,15 +185,15 @@ datasets.
     sampling_times = c(0, 1, 3, 7, 14, 28, 60, 90, 120)
 m_SFO <- mkinmod(parent = mkinsub("SFO"))
 d_SFO_1 <- mkinpredict(m_SFO,
-  c(k_parent_sink = 0.1),
+  c(k_parent = 0.1),
   c(parent = 98), sampling_times)
 d_SFO_1_long <- mkin_wide_to_long(d_SFO_1, time = "time")
 d_SFO_2 <- mkinpredict(m_SFO,
-  c(k_parent_sink = 0.05),
+  c(k_parent = 0.05),
   c(parent = 102), sampling_times)
 d_SFO_2_long <- mkin_wide_to_long(d_SFO_2, time = "time")
 d_SFO_3 <- mkinpredict(m_SFO,
-  c(k_parent_sink = 0.02),
+  c(k_parent = 0.02),
   c(parent = 103), sampling_times)
 d_SFO_3_long <- mkin_wide_to_long(d_SFO_3, time = "time")
 
@@ -214,42 +222,37 @@ datasets.
 #>   Model: value ~ nlme_f(name, time, parent_0, log_k_parent_sink) 
 #>  Data: grouped_data 
 #>        AIC      BIC    logLik
-#>   298.2781 307.7372 -144.1391
+#>   252.7798 262.1358 -121.3899
 #> 
 #> Random effects:
 #>  Formula: list(parent_0 ~ 1, log_k_parent_sink ~ 1)
 #>  Level: ds
 #>  Structure: Diagonal
-#>          parent_0 log_k_parent_sink Residual
-#> StdDev: 0.9374733         0.7098105  3.83543
+#>             parent_0 log_k_parent_sink Residual
+#> StdDev: 0.0006768135         0.6800777 2.489397
 #> 
 #> Fixed effects: parent_0 + log_k_parent_sink ~ 1 
-#>                       Value Std.Error DF  t-value p-value
-#> parent_0          101.76838 1.1445444 45 88.91606       0
-#> log_k_parent_sink  -3.05444 0.4195622 45 -7.28008       0
+#>                       Value Std.Error DF   t-value p-value
+#> parent_0          101.74884 0.6456014 44 157.60321       0
+#> log_k_parent_sink  -3.05575 0.4015811 44  -7.60929       0
 #>  Correlation: 
 #>                   prnt_0
-#> log_k_parent_sink 0.034 
+#> log_k_parent_sink 0.026 
 #> 
 #> Standardized Within-Group Residuals:
 #>        Min         Q1        Med         Q3        Max 
-#> -2.6169360 -0.2185329  0.0574070  0.5720937  3.0459868 
+#> -2.1317488 -0.6878121  0.0828385  0.8592270  2.9529864 
 #> 
-#> Number of Observations: 49
-#> Number of Groups: 3 
plot(augPred(m_nlme, level = 0:1), layout = c(3, 1))
# augPred does not seem to work on fits with more than one state
+#> Number of Observations: 48
+#> Number of Groups: 3 
plot(augPred(m_nlme, level = 0:1), layout = c(3, 1))
# augPred does not seem to work on fits with more than one state
 # variable
 
 
   
-  
-    Contents
-    
-      Arguments
-      Value
-      See also
-      Examples
-    
-
+  
+    
+      Contents
+    
   
 
 
@@ -260,7 +263,7 @@ datasets.
 
 
 
-  Site built with pkgdown 1.4.1.
+  Site built with pkgdown 1.5.1.
 
 
       
-- 
cgit v1.2.1