diff options
| author | Johannes Ranke <jranke@uni-bremen.de> | 2017-11-16 14:55:06 +0100 |
|---|---|---|
| committer | Johannes Ranke <jranke@uni-bremen.de> | 2017-11-16 14:55:06 +0100 |
| commit | 76a0aae725f4d603b3c8e8442bb67081891986b4 (patch) | |
| tree | 1d4f29bd618c00997a9c1247691d9ad45528724e /docs/reference/synthetic_data_for_UBA.html | |
| parent | 545d4bdfe1fc532f6d5d1718ee08444430d2976f (diff) | |
Add the code used for generating synthetic_data_for_UBA
Static documentation except articles rebuilt by pkgdown
Diffstat (limited to 'docs/reference/synthetic_data_for_UBA.html')
| -rw-r--r-- | docs/reference/synthetic_data_for_UBA.html | 290 |
1 files changed, 110 insertions, 180 deletions
diff --git a/docs/reference/synthetic_data_for_UBA.html b/docs/reference/synthetic_data_for_UBA.html index 173af92e..9ff18876 100644 --- a/docs/reference/synthetic_data_for_UBA.html +++ b/docs/reference/synthetic_data_for_UBA.html @@ -70,6 +70,9 @@ <a href="../articles/FOCUS_L.html">Example evaluation of FOCUS Laboratory Data L1 to L3</a> </li> <li> + <a href="../articles/FOCUS_Z.html">Example evaluation of FOCUS Example Dataset Z</a> + </li> + <li> <a href="../articles/compiled_models.html">Performance benefit by using compiled model definitions in mkin</a> </li> <li> @@ -83,12 +86,7 @@ </ul> <ul class="nav navbar-nav navbar-right"> - <li> - <a href="http://github.com/jranke/mkin"> - <span class="fa fa-github fa-lg"></span> - - </a> -</li> + </ul> </div><!--/.nav-collapse --> </div><!--/.container --> @@ -109,7 +107,7 @@ two sequential or two parallel metabolites.</p> <p>Variance component 'a' is based on a normal distribution with standard deviation of 3, Variance component 'b' is also based on a normal distribution, but with a standard deviation of 7. - Variance component 'c' is based on the error model from Rocke and Lorenzato (1995), with the + Variance component 'c' is based on the error model from Rocke and Lorenzato (1995), with the minimum standard deviation (for small y values) of 0.5, and a proportionality constant of 0.07 for the increase of the standard deviation with y.</p> <p>Initial concentrations for metabolites and all values where adding the variance component resulted @@ -133,184 +131,116 @@ <h2 class="hasAnchor" id="source"><a class="anchor" href="#source"></a>Source</h2> <p>Ranke (2014) Prüfung und Validierung von Modellierungssoftware als Alternative - zu ModelMaker 4.0, Umweltbundesamt Projektnummer 27452 - Rocke, David M. und Lorenzato, Stefan (1995) A two-component model for + zu ModelMaker 4.0, Umweltbundesamt Projektnummer 27452</p> +<p>Rocke, David M. und Lorenzato, Stefan (1995) A two-component model for measurement error in analytical chemistry. Technometrics 37(2), 176-184.</p> <h2 class="hasAnchor" id="examples"><a class="anchor" href="#examples"></a>Examples</h2> <pre class="examples"><div class='input'> -<span class='no'>m_synth_SFO_lin</span> <span class='kw'><-</span> <span class='fu'><a href='mkinmod.html'>mkinmod</a></span>(<span class='kw'>parent</span> <span class='kw'>=</span> <span class='fu'>list</span>(<span class='kw'>type</span> <span class='kw'>=</span> <span class='st'>"SFO"</span>, <span class='kw'>to</span> <span class='kw'>=</span> <span class='st'>"M1"</span>), - <span class='kw'>M1</span> <span class='kw'>=</span> <span class='fu'>list</span>(<span class='kw'>type</span> <span class='kw'>=</span> <span class='st'>"SFO"</span>, <span class='kw'>to</span> <span class='kw'>=</span> <span class='st'>"M2"</span>), - <span class='kw'>M2</span> <span class='kw'>=</span> <span class='fu'>list</span>(<span class='kw'>type</span> <span class='kw'>=</span> <span class='st'>"SFO"</span>), <span class='kw'>use_of_ff</span> <span class='kw'>=</span> <span class='st'>"max"</span>)</div><div class='output co'>#> <span class='message'>Successfully compiled differential equation model from auto-generated C code.</span></div><div class='input'> - -<span class='no'>m_synth_SFO_par</span> <span class='kw'><-</span> <span class='fu'><a href='mkinmod.html'>mkinmod</a></span>(<span class='kw'>parent</span> <span class='kw'>=</span> <span class='fu'>list</span>(<span class='kw'>type</span> <span class='kw'>=</span> <span class='st'>"SFO"</span>, <span class='kw'>to</span> <span class='kw'>=</span> <span class='fu'>c</span>(<span class='st'>"M1"</span>, <span class='st'>"M2"</span>), - <span class='kw'>sink</span> <span class='kw'>=</span> <span class='fl'>FALSE</span>), - <span class='kw'>M1</span> <span class='kw'>=</span> <span class='fu'>list</span>(<span class='kw'>type</span> <span class='kw'>=</span> <span class='st'>"SFO"</span>), - <span class='kw'>M2</span> <span class='kw'>=</span> <span class='fu'>list</span>(<span class='kw'>type</span> <span class='kw'>=</span> <span class='st'>"SFO"</span>), <span class='kw'>use_of_ff</span> <span class='kw'>=</span> <span class='st'>"max"</span>)</div><div class='output co'>#> <span class='message'>Successfully compiled differential equation model from auto-generated C code.</span></div><div class='input'> -<span class='no'>m_synth_DFOP_lin</span> <span class='kw'><-</span> <span class='fu'><a href='mkinmod.html'>mkinmod</a></span>(<span class='kw'>parent</span> <span class='kw'>=</span> <span class='fu'>list</span>(<span class='kw'>type</span> <span class='kw'>=</span> <span class='st'>"DFOP"</span>, <span class='kw'>to</span> <span class='kw'>=</span> <span class='st'>"M1"</span>), - <span class='kw'>M1</span> <span class='kw'>=</span> <span class='fu'>list</span>(<span class='kw'>type</span> <span class='kw'>=</span> <span class='st'>"SFO"</span>, <span class='kw'>to</span> <span class='kw'>=</span> <span class='st'>"M2"</span>), - <span class='kw'>M2</span> <span class='kw'>=</span> <span class='fu'>list</span>(<span class='kw'>type</span> <span class='kw'>=</span> <span class='st'>"SFO"</span>), <span class='kw'>use_of_ff</span> <span class='kw'>=</span> <span class='st'>"max"</span>)</div><div class='output co'>#> <span class='message'>Successfully compiled differential equation model from auto-generated C code.</span></div><div class='input'> -<span class='no'>m_synth_DFOP_par</span> <span class='kw'><-</span> <span class='fu'><a href='mkinmod.html'>mkinmod</a></span>(<span class='kw'>parent</span> <span class='kw'>=</span> <span class='fu'>list</span>(<span class='kw'>type</span> <span class='kw'>=</span> <span class='st'>"DFOP"</span>, <span class='kw'>to</span> <span class='kw'>=</span> <span class='fu'>c</span>(<span class='st'>"M1"</span>, <span class='st'>"M2"</span>), - <span class='kw'>sink</span> <span class='kw'>=</span> <span class='fl'>FALSE</span>), - <span class='kw'>M1</span> <span class='kw'>=</span> <span class='fu'>list</span>(<span class='kw'>type</span> <span class='kw'>=</span> <span class='st'>"SFO"</span>), - <span class='kw'>M2</span> <span class='kw'>=</span> <span class='fu'>list</span>(<span class='kw'>type</span> <span class='kw'>=</span> <span class='st'>"SFO"</span>), <span class='kw'>use_of_ff</span> <span class='kw'>=</span> <span class='st'>"max"</span>)</div><div class='output co'>#> <span class='message'>Successfully compiled differential equation model from auto-generated C code.</span></div><div class='input'> -<span class='no'>fit</span> <span class='kw'><-</span> <span class='fu'><a href='mkinfit.html'>mkinfit</a></span>(<span class='no'>m_synth_SFO_lin</span>, <span class='no'>synthetic_data_for_UBA_2014</span><span class='kw'>[[</span><span class='fl'>1</span>]]$<span class='no'>data</span>, <span class='kw'>quiet</span> <span class='kw'>=</span> <span class='fl'>TRUE</span>) -<span class='fu'><a href='plot.mkinfit.html'>plot_sep</a></span>(<span class='no'>fit</span>)</div><img src='synthetic_data_for_UBA_2014-10.png' alt='' width='540' height='400' /><div class='input'><span class='fu'>summary</span>(<span class='no'>fit</span>)</div><div class='output co'>#> mkin version: 0.9.46 -#> R version: 3.4.1 -#> Date of fit: Sat Jul 29 15:15:33 2017 -#> Date of summary: Sat Jul 29 15:15:34 2017 -#> -#> Equations: -#> d_parent/dt = - k_parent * parent -#> d_M1/dt = + f_parent_to_M1 * k_parent * parent - k_M1 * M1 -#> d_M2/dt = + f_M1_to_M2 * k_M1 * M1 - k_M2 * M2 -#> -#> Model predictions using solution type deSolve -#> -#> Fitted with method Port using 381 model solutions performed in 2.241 s -#> -#> Weighting: none -#> -#> Starting values for parameters to be optimised: -#> value type -#> parent_0 101.3500 state -#> k_parent 0.1000 deparm -#> k_M1 0.1001 deparm -#> k_M2 0.1002 deparm -#> f_parent_to_M1 0.5000 deparm -#> f_M1_to_M2 0.5000 deparm -#> -#> Starting values for the transformed parameters actually optimised: -#> value lower upper -#> parent_0 101.350000 -Inf Inf -#> log_k_parent -2.302585 -Inf Inf -#> log_k_M1 -2.301586 -Inf Inf -#> log_k_M2 -2.300587 -Inf Inf -#> f_parent_ilr_1 0.000000 -Inf Inf -#> f_M1_ilr_1 0.000000 -Inf Inf -#> -#> Fixed parameter values: -#> value type -#> M1_0 0 state -#> M2_0 0 state -#> -#> Optimised, transformed parameters with symmetric confidence intervals: -#> Estimate Std. Error Lower Upper -#> parent_0 102.1000 1.71400 98.5800 105.5000 -#> log_k_parent -0.3020 0.04294 -0.3894 -0.2147 -#> log_k_M1 -1.2070 0.07599 -1.3610 -1.0520 -#> log_k_M2 -3.9010 0.06952 -4.0420 -3.7590 -#> f_parent_ilr_1 0.8492 0.18090 0.4812 1.2170 -#> f_M1_ilr_1 0.6780 0.18860 0.2943 1.0620 -#> -#> Parameter correlation: -#> parent_0 log_k_parent log_k_M1 log_k_M2 f_parent_ilr_1 -#> parent_0 1.00000 0.40213 -0.1693 0.02912 -0.4726 -#> log_k_parent 0.40213 1.00000 -0.4210 0.07241 -0.5837 -#> log_k_M1 -0.16931 -0.42103 1.0000 -0.37657 0.7438 -#> log_k_M2 0.02912 0.07241 -0.3766 1.00000 -0.2518 -#> f_parent_ilr_1 -0.47263 -0.58367 0.7438 -0.25177 1.0000 -#> f_M1_ilr_1 0.17148 0.42643 -0.8054 0.52647 -0.8602 -#> f_M1_ilr_1 -#> parent_0 0.1715 -#> log_k_parent 0.4264 -#> log_k_M1 -0.8054 -#> log_k_M2 0.5265 -#> f_parent_ilr_1 -0.8602 -#> f_M1_ilr_1 1.0000 -#> -#> Residual standard error: 2.471 on 33 degrees of freedom -#> -#> 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 102.10000 59.55 1.815e-35 98.58000 105.5000 -#> k_parent 0.73930 23.29 2.337e-22 0.67750 0.8068 -#> k_M1 0.29920 13.16 5.552e-15 0.25630 0.3492 -#> k_M2 0.02023 14.38 4.497e-16 0.01756 0.0233 -#> f_parent_to_M1 0.76870 16.90 4.093e-18 0.66380 0.8483 -#> f_M1_to_M2 0.72290 13.53 2.557e-15 0.60260 0.8178 -#> -#> Chi2 error levels in percent: -#> err.min n.optim df -#> All data 8.454 6 17 -#> parent 8.660 2 6 -#> M1 10.583 2 5 -#> M2 3.586 2 6 -#> -#> Resulting formation fractions: -#> ff -#> parent_M1 0.7687 -#> parent_sink 0.2313 -#> M1_M2 0.7229 -#> M1_sink 0.2771 -#> -#> Estimated disappearance times: -#> DT50 DT90 -#> parent 0.9376 3.114 -#> M1 2.3170 7.697 -#> M2 34.2689 113.839 -#> -#> Data: -#> time variable observed predicted residual -#> 0 parent 101.5 1.021e+02 -0.56248 -#> 0 parent 101.2 1.021e+02 -0.86248 -#> 1 parent 53.9 4.873e+01 5.17118 -#> 1 parent 47.5 4.873e+01 -1.22882 -#> 3 parent 10.4 1.111e+01 -0.70773 -#> 3 parent 7.6 1.111e+01 -3.50773 -#> 7 parent 1.1 5.772e-01 0.52283 -#> 7 parent 0.3 5.772e-01 -0.27717 -#> 14 parent NA 3.264e-03 NA -#> 14 parent 3.5 3.264e-03 3.49674 -#> 28 parent NA 1.045e-07 NA -#> 28 parent 3.2 1.045e-07 3.20000 -#> 60 parent NA -1.054e-10 NA -#> 60 parent NA -1.054e-10 NA -#> 90 parent 0.6 -1.875e-11 0.60000 -#> 90 parent NA -1.875e-11 NA -#> 120 parent NA -2.805e-11 NA -#> 120 parent 3.5 -2.805e-11 3.50000 -#> 0 M1 NA 0.000e+00 NA -#> 0 M1 NA 0.000e+00 NA -#> 1 M1 36.4 3.479e+01 1.61088 -#> 1 M1 37.4 3.479e+01 2.61088 -#> 3 M1 34.3 3.937e+01 -5.07027 -#> 3 M1 39.8 3.937e+01 0.42973 -#> 7 M1 15.1 1.549e+01 -0.38715 -#> 7 M1 17.8 1.549e+01 2.31285 -#> 14 M1 5.8 1.995e+00 3.80469 -#> 14 M1 1.2 1.995e+00 -0.79531 -#> 28 M1 NA 3.034e-02 NA -#> 28 M1 NA 3.034e-02 NA -#> 60 M1 0.5 2.111e-06 0.50000 -#> 60 M1 NA 2.111e-06 NA -#> 90 M1 NA 2.913e-10 NA -#> 90 M1 3.2 2.913e-10 3.20000 -#> 120 M1 1.5 3.625e-11 1.50000 -#> 120 M1 0.6 3.625e-11 0.60000 -#> 0 M2 NA 0.000e+00 NA -#> 0 M2 NA 0.000e+00 NA -#> 1 M2 NA 4.455e+00 NA -#> 1 M2 4.8 4.455e+00 0.34517 -#> 3 M2 20.9 2.153e+01 -0.62527 -#> 3 M2 19.3 2.153e+01 -2.22527 -#> 7 M2 42.0 4.192e+01 0.07941 -#> 7 M2 43.1 4.192e+01 1.17941 -#> 14 M2 49.4 4.557e+01 3.83353 -#> 14 M2 44.3 4.557e+01 -1.26647 -#> 28 M2 34.6 3.547e+01 -0.87275 -#> 28 M2 33.0 3.547e+01 -2.47275 -#> 60 M2 18.8 1.858e+01 0.21837 -#> 60 M2 17.6 1.858e+01 -0.98163 -#> 90 M2 10.6 1.013e+01 0.47130 -#> 90 M2 10.8 1.013e+01 0.67130 -#> 120 M2 9.8 5.521e+00 4.27893 -#> 120 M2 3.3 5.521e+00 -2.22107</div><div class='input'> -</div></pre> +# The data have been generated using the following kinetic models +m_synth_SFO_lin <- mkinmod(parent = list(type = "SFO", to = "M1"), + M1 = list(type = "SFO", to = "M2"), + M2 = list(type = "SFO"), use_of_ff = "max") + + +m_synth_SFO_par <- mkinmod(parent = list(type = "SFO", to = c("M1", "M2"), + sink = FALSE), + M1 = list(type = "SFO"), + M2 = list(type = "SFO"), use_of_ff = "max") + +m_synth_DFOP_lin <- mkinmod(parent = list(type = "DFOP", to = "M1"), + M1 = list(type = "SFO", to = "M2"), + M2 = list(type = "SFO"), use_of_ff = "max") + +m_synth_DFOP_par <- mkinmod(parent = list(type = "DFOP", to = c("M1", "M2"), + sink = FALSE), + M1 = list(type = "SFO"), + M2 = list(type = "SFO"), use_of_ff = "max") + +# The model predictions without intentional error were generated as follows +sampling_times = c(0, 1, 3, 7, 14, 28, 60, 90, 120) + +d_synth_SFO_lin <- mkinpredict(m_synth_SFO_lin, + c(k_parent = 0.7, f_parent_to_M1 = 0.8, + k_M1 = 0.3, f_M1_to_M2 = 0.7, + k_M2 = 0.02), + c(parent = 100, M1 = 0, M2 = 0), + sampling_times) + +d_synth_DFOP_lin <- mkinpredict(m_synth_DFOP_lin, + c(k1 = 0.2, k2 = 0.02, g = 0.5, + f_parent_to_M1 = 0.5, k_M1 = 0.3, + f_M1_to_M2 = 0.7, k_M2 = 0.02), + c(parent = 100, M1 = 0, M2 = 0), + sampling_times) + +d_synth_SFO_par <- mkinpredict(m_synth_SFO_par, + c(k_parent = 0.2, + f_parent_to_M1 = 0.8, k_M1 = 0.01, + f_parent_to_M2 = 0.2, k_M2 = 0.02), + c(parent = 100, M1 = 0, M2 = 0), + sampling_times) + +d_synth_DFOP_par <- mkinpredict(m_synth_DFOP_par, + c(k1 = 0.3, k2 = 0.02, g = 0.7, + f_parent_to_M1 = 0.6, k_M1 = 0.04, + f_parent_to_M2 = 0.4, k_M2 = 0.01), + c(parent = 100, M1 = 0, M2 = 0), + sampling_times) + +# Construct names for datasets with errors +d_synth_names = paste0("d_synth_", c("SFO_lin", "SFO_par", + "DFOP_lin", "DFOP_par")) + +# Function for adding errors. The add_err function now published with this +# package is a slightly generalised version where the names of secondary +# compartments that should have an initial value of zero (M1 and M2 in this +# case) are not hardcoded any more. +add_err = function(d, sdfunc, LOD = 0.1, reps = 2, seed = 123456789) +{ + set.seed(seed) + d_long = mkin_wide_to_long(d, time = "time") + d_rep = data.frame(lapply(d_long, rep, each = 2)) + d_rep$value = rnorm(length(d_rep$value), d_rep$value, sdfunc(d_rep$value)) + + d_rep[d_rep$time == 0 & d_rep$name <!-- %in% c("M1", "M2"), "value"] <- 0 --> + d_NA <- transform(d_rep, value = ifelse(value < LOD, NA, value)) + d_NA$value <- round(d_NA$value, 1) + return(d_NA) +} + +# The following is the two-component model of Rocke and Lorenzato (1995) +sdfunc_twocomp = function(value, sd_low, rsd_high) { + sqrt(sd_low^2 + value^2 * rsd_high^2) +} + +# Add the errors. +for (d_synth_name in d_synth_names) +{ + d_synth = get(d_synth_name) + assign(paste0(d_synth_name, "_a"), add_err(d_synth, function(value) 3)) + assign(paste0(d_synth_name, "_b"), add_err(d_synth, function(value) 7)) + assign(paste0(d_synth_name, "_c"), add_err(d_synth, + function(value) sdfunc_twocomp(value, 0.5, 0.07))) + +} + +d_synth_err_names = c( + paste(rep(d_synth_names, each = 3), letters[1:3], sep = "_") +) + +# This is just one example of an evaluation using the kinetic model used for +# the generation of the data +fit <- mkinfit(m_synth_SFO_lin, synthetic_data_for_UBA_2014[[1]]$data, + quiet = TRUE) +plot_sep(fit) +summary(fit) + +</div><div class='output co'>#> <span class='error'>Error: <text>:68:43: Unerwartete(s) SPECIAL</span> +#> <span class='error'>67: </span> +#> <span class='error'>68: d_rep[d_rep$time == 0 & d_rep$name <!-- %in%</span> +#> <span class='error'> ^</span></div></pre> </div> <div class="col-md-3 hidden-xs hidden-sm" id="sidebar"> <h2>Contents</h2> |