From 6733555d7a9315c55001770bacc4c61c4d4f39d5 Mon Sep 17 00:00:00 2001 From: Johannes Ranke Date: Sun, 21 Jun 2015 01:46:51 +0200 Subject: Do the t-test for untransformed parameters --- .gitignore | 3 +- DESCRIPTION | 4 +- NEWS.md | 1 + R/mkinfit.R | 89 +++++-- README.Rmd | 2 + README.md | 2 + figure/unnamed-chunk-7-2.png | Bin 6533 -> 6643 bytes tests/testthat/test_FOCUS_D_UBA_expertise.R | 28 +++ vignettes/FOCUS_D.html | 50 ++-- vignettes/FOCUS_L.html | 357 +++++++++++++++------------- vignettes/FOCUS_Z.Rnw | 2 + vignettes/FOCUS_Z.pdf | Bin 215014 -> 224256 bytes vignettes/compiled_models.html | 42 ++-- vignettes/mkin.pdf | Bin 160260 -> 160260 bytes 14 files changed, 340 insertions(+), 240 deletions(-) diff --git a/.gitignore b/.gitignore index fe72b773..a06689ac 100644 --- a/.gitignore +++ b/.gitignore @@ -11,7 +11,8 @@ vignettes/*.R vignettes/.build.timestamp vignettes/mkin.tex vignettes/FOCUS_Z.tex +vignettes/FOCUS_Z.fls +vignettes/FOCUS_Z.fdb_latexmk vignettes/*cache/ -======= vignettes/compiled_models_cache mkin.Rcheck diff --git a/DESCRIPTION b/DESCRIPTION index b5045e4e..6586dd8b 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -3,7 +3,7 @@ Type: Package Title: Routines for Fitting Kinetic Models with One or More State Variables to Chemical Degradation Data Version: 0.9-36 -Date: 2015-06-20 +Date: 2015-06-21 Authors@R: c(person("Johannes", "Ranke", role = c("aut", "cre", "cph"), email = "jranke@uni-bremen.de"), person("Katrin", "Lindenberger", role = "ctb"), @@ -13,7 +13,7 @@ Description: Calculation routines based on the FOCUS Kinetics Report (2006). Includes a function for conveniently defining differential equation models, model solution based on eigenvalues if possible or using numerical solvers and a choice of the optimisation methods made available by the 'FME' package. - If a C compiler (on windows: 'Rtools') are installed, differential + If a C compiler (on windows: 'Rtools') is installed, differential equation models are solved using compiled C functions. Please note that no warranty is implied for correctness of results or fitness for a particular purpose. diff --git a/NEWS.md b/NEWS.md index 1256c652..7ee9696f 100644 --- a/NEWS.md +++ b/NEWS.md @@ -2,6 +2,7 @@ ## MAJOR CHANGES +- `summary.mkinfit()`: A one-sided t-test for significant difference of untransformed parameters from zero is now always shown, based on the assumption of normal distribution for estimators of all untransformed parameters. Use with caution, as this assumption is unrealistic e.g. for rate constants in these nonlinear kinetic models. - If a compiler (gcc) is installed, use a version of the differential equation model compiled from C code, which is a huge performance boost for models where only the deSolve method works. - `mkinmod()`: Create a list component $cf (of class CFuncList) in the list returned by mkinmod, if a version can be compiled from autogenerated C code (see above). - `mkinfit()`: Set the default `solution_type` to `deSolve` when a compiled version of the model is present, except when an analytical solution is possible. diff --git a/R/mkinfit.R b/R/mkinfit.R index 7f9a55d9..cfd87829 100644 --- a/R/mkinfit.R +++ b/R/mkinfit.R @@ -234,10 +234,19 @@ mkinfit <- function(mkinmod, observed, } } + # Define outtimes for model solution. + # Include time points at which observed data are available + # Make sure we include time 0, so initial values for state variables are for time 0 + outtimes = sort(unique(c(observed$time, seq(min(observed$time), + max(observed$time), + length.out = n.outtimes)))) + + cost.old <- 1e100 # The first model cost should be smaller than this value calls <- 0 # Counter for number of model solutions out_predicted <- NA - # Define the model cost function + + # Define the model cost function for optimisation, including (back)transformations cost <- function(P) { assign("calls", calls+1, inherits=TRUE) # Increase the model solution counter @@ -245,12 +254,6 @@ mkinfit <- function(mkinmod, observed, # Trace parameter values if requested if(trace_parms) cat(P, "\n") - # Time points at which observed data are available - # Make sure we include time 0, so initial values for state variables are for time 0 - outtimes = sort(unique(c(observed$time, seq(min(observed$time), - max(observed$time), - length.out = n.outtimes)))) - if(length(state.ini.optim) > 0) { odeini <- c(P[1:length(state.ini.optim)], state.ini.fixed) names(odeini) <- c(state.ini.optim.boxnames, state.ini.fixed.boxnames) @@ -313,6 +316,35 @@ mkinfit <- function(mkinmod, observed, return(mC) } + # Define the model cost function for the t-test, without parameter transformation + cost_notrans <- function(P) + { + if(length(state.ini.optim) > 0) { + odeini <- c(P[1:length(state.ini.optim)], state.ini.fixed) + names(odeini) <- c(state.ini.optim.boxnames, state.ini.fixed.boxnames) + } else { + odeini <- state.ini.fixed + names(odeini) <- state.ini.fixed.boxnames + } + + odeparms <- c(P[(length(state.ini.optim) + 1):length(P)], parms.fixed) + + # Solve the system with current transformed parameter values + out <- mkinpredict(mkinmod, odeparms, + odeini, outtimes, + solution_type = solution_type, + use_compiled = use_compiled, + method.ode = method.ode, + atol = atol, rtol = rtol, ...) + + mC <- modCost(out, observed, y = "value", + err = err, weight = weight, scaleVar = scaleVar) + + return(mC) + } + + # Define lower and upper bounds other than -Inf and Inf for parameters + # for which no internal transformation is requested in the call to mkinfit. lower <- rep(-Inf, length(c(state.ini.optim, transparms.optim))) upper <- rep(Inf, length(c(state.ini.optim, transparms.optim))) names(lower) <- names(upper) <- c(names(state.ini.optim), names(transparms.optim)) @@ -431,6 +463,17 @@ mkinfit <- function(mkinmod, observed, bparms.fixed = c(state.ini.fixed, parms.fixed) bparms.all = c(bparms.optim, parms.fixed) + # Attach the cost functions to the fit for post-hoc parameter uncertainty analysis + fit$cost <- cost + fit$cost_notrans <- cost_notrans + + # Estimate the Hessian for the model cost without parameter transformations + # to make it possible to obtain the usual t-test + # Code ported from FME::modFit + Jac_notrans <- gradient(function(p, ...) cost_notrans(p)$residuals$res, + bparms.optim, centered = TRUE) + fit$hessian_notrans <- 2 * t(Jac_notrans) %*% Jac_notrans + # Collect initial parameter values in three dataframes fit$start <- data.frame(value = c(state.ini.optim, parms.optim)) @@ -489,29 +532,36 @@ summary.mkinfit <- function(object, data = TRUE, distimes = TRUE, alpha = 0.05, p <- length(param) mod_vars <- names(object$mkinmod$diffs) covar <- try(solve(0.5*object$hessian), silent = TRUE) # unscaled covariance + covar_notrans <- try(solve(0.5*object$hessian_notrans), silent = TRUE) # unscaled covariance rdf <- object$df.residual resvar <- object$ssr / rdf if (!is.numeric(covar)) { covar <- NULL - se <- lci <- uci <- tval <- pval1 <- pval2 <- rep(NA, p) + se <- lci <- uci <- rep(NA, p) } else { rownames(covar) <- colnames(covar) <- pnames se <- sqrt(diag(covar) * resvar) lci <- param + qt(alpha/2, rdf) * se uci <- param + qt(1-alpha/2, rdf) * se - tval <- param/se - pval1 <- 2 * pt(abs(tval), rdf, lower.tail = FALSE) - pval2 <- pt(abs(tval), rdf, lower.tail = FALSE) + } + + if (!is.numeric(covar_notrans)) { + covar_notrans <- NULL + se_notrans <- tval <- pval <- rep(NA, p) + } else { + rownames(covar_notrans) <- colnames(covar_notrans) <- bpnames + se_notrans <- sqrt(diag(covar_notrans) * resvar) + tval <- object$bparms.optim/se_notrans + pval <- pt(abs(tval), rdf, lower.tail = FALSE) } names(se) <- pnames modVariance <- object$ssr / length(object$residuals) - param <- cbind(param, se, lci, uci, tval, pval1, pval2) - dimnames(param) <- list(pnames, c("Estimate", "Std. Error", "Lower", "Upper", - "t value", "Pr(>|t|)", "Pr(>t)")) + param <- cbind(param, se, lci, uci) + dimnames(param) <- list(pnames, c("Estimate", "Std. Error", "Lower", "Upper")) - bparam <- cbind(Estimate = object$bparms.optim, Lower = NA, Upper = NA) + bparam <- cbind(Estimate = object$bparms.optim, se_notrans, "t value" = tval, "Pr(>t)" = pval, Lower = NA, Upper = NA) # Transform boundaries of CI for one parameter at a time, # with the exception of sets of formation fractions (single fractions are OK). @@ -617,7 +667,7 @@ print.summary.mkinfit <- function(x, digits = max(3, getOption("digits") - 3), . if(length(x$fixed$value) == 0) cat("None\n") else print(x$fixed) - cat("\nOptimised, transformed parameters:\n") + cat("\nOptimised, transformed parameters with symmetric confidence intervals:\n") print(signif(x$par, digits = digits)) if (x$calls > 0) { @@ -634,8 +684,11 @@ print.summary.mkinfit <- function(x, digits = max(3, getOption("digits") - 3), . cat("\nResidual standard error:", format(signif(x$sigma, digits)), "on", rdf, "degrees of freedom\n") - cat("\nBacktransformed parameters:\n") - print(signif(x$bpar, digits = digits)) + cat("\nBacktransformed parameters:\n", + " Confidence intervals for internally transformed parameters are asymmetric.\n", + " t-test (unrealistically) based on the assumption of normal distribution\n", + " for estimators of untransformed parameters.\n", sep = "") + print(signif(x$bpar[, c(1, 3, 4, 5, 6)], digits = digits)) cat("\nChi2 error levels in percent:\n") x$errmin$err.min <- 100 * x$errmin$err.min diff --git a/README.Rmd b/README.Rmd index c5df2aef..e7cc5798 100644 --- a/README.Rmd +++ b/README.Rmd @@ -130,6 +130,8 @@ documentation or the package vignettes referenced from the not include meaningless values like negative rate constants or formation fractions adding up to more than 1, which can not occur in a single experiment with a single defined radiolabel position. +* The usual one-sided t-test for significant difference from zero is nevertheless + shown based on estimators for the untransformed parameters. * Summary and plotting functions. The `summary` of an `mkinfit` object is in fact a full report that should give enough information to be able to approximately reproduce the fit with other tools. diff --git a/README.md b/README.md index 37af4f31..92f1e8b1 100644 --- a/README.md +++ b/README.md @@ -144,6 +144,8 @@ documentation or the package vignettes referenced from the not include meaningless values like negative rate constants or formation fractions adding up to more than 1, which can not occur in a single experiment with a single defined radiolabel position. +* The usual one-sided t-test for significant difference from zero is nevertheless + shown based on estimators for the untransformed parameters. * Summary and plotting functions. The `summary` of an `mkinfit` object is in fact a full report that should give enough information to be able to approximately reproduce the fit with other tools. diff --git a/figure/unnamed-chunk-7-2.png b/figure/unnamed-chunk-7-2.png index b4058df9..c31a5f53 100644 Binary files a/figure/unnamed-chunk-7-2.png and b/figure/unnamed-chunk-7-2.png differ diff --git a/tests/testthat/test_FOCUS_D_UBA_expertise.R b/tests/testthat/test_FOCUS_D_UBA_expertise.R index b12834e5..5d5a801a 100644 --- a/tests/testthat/test_FOCUS_D_UBA_expertise.R +++ b/tests/testthat/test_FOCUS_D_UBA_expertise.R @@ -56,3 +56,31 @@ test_that("DT50/90 are correct for FOCUS D when not using formation fractions", # References: # Ranke (2014) Prüfung und Validierung von Modellierungssoftware als Alternative # zu ModelMaker 4.0, Umweltbundesamt Projektnummer 27452 + +context("The t-test for significant difference from zero") + +test_that("The t-value for fits using internal transformations corresponds with result from FME", { + + expect_equal(signif(summary(fit.default)$bpar[, "t value"], 5), + c(parent_0 = 61.720, k_parent_sink = 12.777, k_parent_m1 = 24.248, k_m1_sink = 7.3486)) + +}) + +m_synth_DFOP_par.minff <- mkinmod(parent = mkinsub("DFOP", c("M1", "M2")), + M1 = mkinsub("SFO"), + M2 = mkinsub("SFO"), + use_of_ff = "min", quiet = TRUE) + +fit_DFOP_par_c_2 <- mkinfit(m_synth_DFOP_par.minff, + synthetic_data_for_UBA_2014[[12]]$data, + quiet = TRUE) + +test_that("The t-value for fits using internal transformations corresponds with results from FME, synthetic data", { + + # Note that the k1 and k2 are exchanged in the untransformed fit evaluated with FME for this test + expect_equal(signif(summary(fit_DFOP_par_c_2)$bpar[1:7, "t value"], 5), + c(parent_0 = 80.054, k_M1_sink = 12.291, k_M2_sink = 10.588, + f_parent_to_M1 = 21.4960, f_parent_to_M2 = 24.0890, + k1 = 16.1450, k2 = 8.1747)) + +}) diff --git a/vignettes/FOCUS_D.html b/vignettes/FOCUS_D.html index 6573cc7a..b1ea64ea 100644 --- a/vignettes/FOCUS_D.html +++ b/vignettes/FOCUS_D.html @@ -215,13 +215,7 @@ library we look a the data. We have observed concentrations in the column named named parent and m1.

library("mkin")
-
- -
## Loading required package: minpack.lm
-## Loading required package: rootSolve
-
- -
print(FOCUS_2006_D)
+print(FOCUS_2006_D)
 
##      name time  value
@@ -276,7 +270,7 @@ kinetics (SFO) to one metabolite named m1, which also degrades with SFO kinetics
 
 

The call to mkinmod returns a degradation model. The differential equations represented in R code can be found in the character vector $diffs of the mkinmod object. If -the ccSolve package is installed and functional, the differential equation model +the gcc compiler is installed and functional, the differential equation model will be compiled from auto-generated C code.

SFO_SFO <- mkinmod(parent = mkinsub("SFO", "m1"), m1 = mkinsub("SFO"))
@@ -312,7 +306,7 @@ using the plot method for mkinfit objects.

mkinparplot(fit)
 
-

plot of chunk unnamed-chunk-6

+

plot of chunk unnamed-chunk-6

A comprehensive report of the results is obtained using the summary method for mkinfit objects.

@@ -321,9 +315,9 @@ objects.

## mkin version:    0.9.36 
-## R version:       3.2.0 
-## Date of fit:     Fri Jun  5 14:20:31 2015 
-## Date of summary: Fri Jun  5 14:20:31 2015 
+## R version:       3.2.1 
+## Date of fit:     Sun Jun 21 01:47:59 2015 
+## Date of summary: Sun Jun 21 01:47:59 2015 
 ## 
 ## Equations:
 ## d_parent = - k_parent_sink * parent - k_parent_m1 * parent
@@ -331,7 +325,7 @@ objects.

## ## Model predictions using solution type deSolve ## -## Fitted with method Port using 153 model solutions performed in 0.621 s +## Fitted with method Port using 153 model solutions performed in 0.698 s ## ## Weighting: none ## @@ -353,17 +347,12 @@ objects.

## value type ## m1_0 0 state ## -## Optimised, transformed parameters: -## Estimate Std. Error Lower Upper t value Pr(>|t|) -## parent_0 99.600 1.61400 96.330 102.900 61.72 4.048e-38 -## log_k_parent_sink -3.038 0.07826 -3.197 -2.879 -38.82 5.601e-31 -## log_k_parent_m1 -2.980 0.04124 -3.064 -2.897 -72.27 1.446e-40 -## log_k_m1_sink -5.248 0.13610 -5.523 -4.972 -38.56 7.087e-31 -## Pr(>t) -## parent_0 2.024e-38 -## log_k_parent_sink 2.800e-31 -## log_k_parent_m1 7.228e-41 -## log_k_m1_sink 3.543e-31 +## Optimised, transformed parameters with symmetric confidence intervals: +## Estimate Std. Error Lower Upper +## parent_0 99.600 1.61400 96.330 102.900 +## log_k_parent_sink -3.038 0.07826 -3.197 -2.879 +## log_k_parent_m1 -2.980 0.04124 -3.064 -2.897 +## log_k_m1_sink -5.248 0.13610 -5.523 -4.972 ## ## Parameter correlation: ## parent_0 log_k_parent_sink log_k_parent_m1 log_k_m1_sink @@ -375,11 +364,14 @@ objects.

## Residual standard error: 3.211 on 36 degrees of freedom ## ## Backtransformed parameters: -## Estimate Lower Upper -## parent_0 99.600000 96.330000 1.029e+02 -## k_parent_sink 0.047920 0.040890 5.616e-02 -## k_parent_m1 0.050780 0.046700 5.521e-02 -## k_m1_sink 0.005261 0.003992 6.933e-03 +## 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 99.600000 61.720 2.024e-38 96.330000 1.029e+02 +## k_parent_sink 0.047920 12.780 3.050e-15 0.040890 5.616e-02 +## k_parent_m1 0.050780 24.250 3.407e-24 0.046700 5.521e-02 +## k_m1_sink 0.005261 7.349 5.758e-09 0.003992 6.933e-03 ## ## Chi2 error levels in percent: ## err.min n.optim df diff --git a/vignettes/FOCUS_L.html b/vignettes/FOCUS_L.html index 96ea70ce..692caf93 100644 --- a/vignettes/FOCUS_L.html +++ b/vignettes/FOCUS_L.html @@ -214,13 +214,7 @@ hr { report, p. 284:

library("mkin")
-
- -
## Loading required package: minpack.lm
-## Loading required package: rootSolve
-
- -
FOCUS_2006_L1 = data.frame(
+FOCUS_2006_L1 = data.frame(
   t = rep(c(0, 1, 2, 3, 5, 7, 14, 21, 30), each = 2),
   parent = c(88.3, 91.4, 85.6, 84.5, 78.9, 77.6, 
              72.0, 71.9, 50.3, 59.4, 47.0, 45.1,
@@ -242,17 +236,17 @@ given in the FOCUS report. 

summary(m.L1.SFO)
-
## mkin version:    0.9.35 
-## R version:       3.1.2 
-## Date of fit:     Sat Feb 21 14:44:53 2015 
-## Date of summary: Sat Feb 21 14:44:53 2015 
+
## mkin version:    0.9.36 
+## R version:       3.2.1 
+## Date of fit:     Sun Jun 21 01:47:59 2015 
+## Date of summary: Sun Jun 21 01:47:59 2015 
 ## 
 ## Equations:
 ## d_parent = - k_parent_sink * parent
 ## 
 ## Model predictions using solution type analytical 
 ## 
-## Fitted with method Port using 37 model solutions performed in 0.098 s
+## Fitted with method Port using 37 model solutions performed in 0.093 s
 ## 
 ## Weighting: none
 ## 
@@ -269,13 +263,10 @@ summary(m.L1.SFO)
 ## Fixed parameter values:
 ## None
 ## 
-## Optimised, transformed parameters:
-##                   Estimate Std. Error  Lower  Upper t value  Pr(>|t|)
-## parent_0            92.470    1.36800 89.570 95.370   67.58 4.339e-21
-## log_k_parent_sink   -2.347    0.04057 -2.433 -2.261  -57.86 5.155e-20
-##                      Pr(>t)
-## parent_0          2.170e-21
-## log_k_parent_sink 2.577e-20
+## Optimised, transformed parameters with symmetric confidence intervals:
+##                   Estimate Std. Error  Lower  Upper
+## parent_0            92.470    1.36800 89.570 95.370
+## log_k_parent_sink   -2.347    0.04057 -2.433 -2.261
 ## 
 ## Parameter correlation:
 ##                   parent_0 log_k_parent_sink
@@ -285,9 +276,12 @@ summary(m.L1.SFO)
 ## Residual standard error: 2.948 on 16 degrees of freedom
 ## 
 ## Backtransformed parameters:
-##               Estimate    Lower   Upper
-## parent_0      92.47000 89.57000 95.3700
-## k_parent_sink  0.09561  0.08773  0.1042
+##   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      92.47000   67.58 2.170e-21 89.57000 95.3700
+## k_parent_sink  0.09561   24.65 1.867e-14  0.08773  0.1042
 ## 
 ## Chi2 error levels in percent:
 ##          err.min n.optim df
@@ -341,20 +335,31 @@ The residual plot can be easily obtained by

is checked.

m.L1.FOMC <- mkinfit("FOMC", FOCUS_2006_L1_mkin, quiet=TRUE)
-summary(m.L1.FOMC, data = FALSE)
 
-
## mkin version:    0.9.35 
-## R version:       3.1.2 
-## Date of fit:     Sat Feb 21 14:44:55 2015 
-## Date of summary: Sat Feb 21 14:44:55 2015 
+
## Warning in mkinfit("FOMC", FOCUS_2006_L1_mkin, quiet = TRUE): Optimisation by method Port did not converge.
+## Convergence code is 1
+
+ +
summary(m.L1.FOMC, data = FALSE)
+
+ +
## mkin version:    0.9.36 
+## R version:       3.2.1 
+## Date of fit:     Sun Jun 21 01:48:00 2015 
+## Date of summary: Sun Jun 21 01:48:00 2015 
+## 
+## 
+## Warning: Optimisation by method Port did not converge.
+## Convergence code is 1 
+## 
 ## 
 ## Equations:
-## d_parent = - (alpha/beta) * ((time/beta) + 1)^-1 * parent
+## d_parent = - (alpha/beta) * 1/((time/beta) + 1) * parent
 ## 
 ## Model predictions using solution type analytical 
 ## 
-## Fitted with method Port using 611 model solutions performed in 1.509 s
+## Fitted with method Port using 188 model solutions performed in 0.463 s
 ## 
 ## Weighting: none
 ## 
@@ -373,29 +378,28 @@ summary(m.L1.FOMC, data = FALSE)
 ## Fixed parameter values:
 ## None
 ## 
-## Optimised, transformed parameters:
-##           Estimate Std. Error    Lower   Upper  t value  Pr(>|t|)
-## parent_0     92.47      1.482    89.31   95.63 62.39000 1.546e-19
-## log_alpha    11.25    598.200 -1264.00 1286.00  0.01880 9.852e-01
-## log_beta     13.60    598.200 -1261.00 1289.00  0.02273 9.822e-01
-##              Pr(>t)
-## parent_0  7.730e-20
-## log_alpha 4.926e-01
-## log_beta  4.911e-01
+## Optimised, transformed parameters with symmetric confidence intervals:
+##           Estimate Std. Error  Lower Upper
+## parent_0     92.47      1.422  89.44 95.50
+## log_alpha    15.43     15.080 -16.71 47.58
+## log_beta     17.78     15.090 -14.37 49.93
 ## 
 ## Parameter correlation:
 ##           parent_0 log_alpha log_beta
-## parent_0    1.0000   -0.3016  -0.3016
-## log_alpha  -0.3016    1.0000   1.0000
-## log_beta   -0.3016    1.0000   1.0000
+## parent_0    1.0000    0.1129   0.1112
+## log_alpha   0.1129    1.0000   1.0000
+## log_beta    0.1112    1.0000   1.0000
 ## 
 ## Residual standard error: 3.045 on 15 degrees of freedom
 ## 
 ## Backtransformed parameters:
-##           Estimate Lower Upper
-## parent_0     92.47 89.31 95.63
-## alpha     76830.00  0.00   Inf
-## beta     803500.00  0.00   Inf
+##   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 9.247e+01  65.150 4.044e-20 8.944e+01 9.550e+01
+## alpha    5.044e+06   1.271 1.115e-01 5.510e-08 4.618e+20
+## beta     5.276e+07   1.259 1.137e-01 5.732e-07 4.857e+21
 ## 
 ## Chi2 error levels in percent:
 ##          err.min n.optim df
@@ -403,8 +407,8 @@ summary(m.L1.FOMC, data = FALSE)
 ## parent     3.619       3  6
 ## 
 ## Estimated disappearance times:
-##         DT50  DT90 DT50back
-## parent 7.249 24.08    7.249
+##        DT50  DT90 DT50back
+## parent 7.25 24.08     7.25
 

Due to the higher number of parameters, and the lower number of degrees of @@ -442,17 +446,17 @@ FOCUS_2006_L2_mkin <- mkin_wide_to_long(FOCUS_2006_L2) summary(m.L2.SFO)

-
## mkin version:    0.9.35 
-## R version:       3.1.2 
-## Date of fit:     Sat Feb 21 14:44:55 2015 
-## Date of summary: Sat Feb 21 14:44:55 2015 
+
## mkin version:    0.9.36 
+## R version:       3.2.1 
+## Date of fit:     Sun Jun 21 01:48:00 2015 
+## Date of summary: Sun Jun 21 01:48:00 2015 
 ## 
 ## Equations:
 ## d_parent = - k_parent_sink * parent
 ## 
 ## Model predictions using solution type analytical 
 ## 
-## Fitted with method Port using 41 model solutions performed in 0.1 s
+## Fitted with method Port using 41 model solutions performed in 0.097 s
 ## 
 ## Weighting: none
 ## 
@@ -469,13 +473,10 @@ summary(m.L2.SFO)
 ## Fixed parameter values:
 ## None
 ## 
-## Optimised, transformed parameters:
-##                   Estimate Std. Error   Lower   Upper t value  Pr(>|t|)
-## parent_0           91.4700     3.8070 82.9800 99.9500  24.030 3.545e-10
-## log_k_parent_sink  -0.4112     0.1074 -0.6505 -0.1719  -3.828 3.329e-03
-##                      Pr(>t)
-## parent_0          1.773e-10
-## log_k_parent_sink 1.664e-03
+## Optimised, transformed parameters with symmetric confidence intervals:
+##                   Estimate Std. Error   Lower   Upper
+## parent_0           91.4700     3.8070 82.9800 99.9500
+## log_k_parent_sink  -0.4112     0.1074 -0.6505 -0.1719
 ## 
 ## Parameter correlation:
 ##                   parent_0 log_k_parent_sink
@@ -485,9 +486,12 @@ summary(m.L2.SFO)
 ## Residual standard error: 5.51 on 10 degrees of freedom
 ## 
 ## Backtransformed parameters:
-##               Estimate   Lower   Upper
-## parent_0       91.4700 82.9800 99.9500
-## k_parent_sink   0.6629  0.5218  0.8421
+##   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       91.4700   24.03 1.773e-10 82.9800 99.9500
+## k_parent_sink   0.6629    9.31 1.525e-06  0.5218  0.8421
 ## 
 ## Chi2 error levels in percent:
 ##          err.min n.optim df
@@ -552,17 +556,17 @@ mkinresplot(m.L2.FOMC)
 
summary(m.L2.FOMC, data = FALSE)
 
-
## mkin version:    0.9.35 
-## R version:       3.1.2 
-## Date of fit:     Sat Feb 21 14:44:55 2015 
-## Date of summary: Sat Feb 21 14:44:55 2015 
+
## mkin version:    0.9.36 
+## R version:       3.2.1 
+## Date of fit:     Sun Jun 21 01:48:00 2015 
+## Date of summary: Sun Jun 21 01:48:00 2015 
 ## 
 ## Equations:
-## d_parent = - (alpha/beta) * ((time/beta) + 1)^-1 * parent
+## d_parent = - (alpha/beta) * 1/((time/beta) + 1) * parent
 ## 
 ## Model predictions using solution type analytical 
 ## 
-## Fitted with method Port using 81 model solutions performed in 0.201 s
+## Fitted with method Port using 81 model solutions performed in 0.191 s
 ## 
 ## Weighting: none
 ## 
@@ -581,11 +585,11 @@ mkinresplot(m.L2.FOMC)
 ## Fixed parameter values:
 ## None
 ## 
-## Optimised, transformed parameters:
-##           Estimate Std. Error   Lower   Upper t value  Pr(>|t|)    Pr(>t)
-## parent_0   93.7700     1.8560 89.5700 97.9700 50.5100 2.345e-12 1.173e-12
-## log_alpha   0.3180     0.1867 -0.1044  0.7405  1.7030 1.227e-01 6.137e-02
-## log_beta    0.2102     0.2943 -0.4555  0.8759  0.7142 4.932e-01 2.466e-01
+## Optimised, transformed parameters with symmetric confidence intervals:
+##           Estimate Std. Error   Lower   Upper
+## parent_0   93.7700     1.8560 89.5700 97.9700
+## log_alpha   0.3180     0.1867 -0.1044  0.7405
+## log_beta    0.2102     0.2943 -0.4555  0.8759
 ## 
 ## Parameter correlation:
 ##           parent_0 log_alpha log_beta
@@ -596,10 +600,13 @@ mkinresplot(m.L2.FOMC)
 ## Residual standard error: 2.628 on 9 degrees of freedom
 ## 
 ## Backtransformed parameters:
-##          Estimate   Lower  Upper
-## parent_0   93.770 89.5700 97.970
-## alpha       1.374  0.9009  2.097
-## beta        1.234  0.6341  2.401
+##   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   93.770  50.510 1.173e-12 89.5700 97.970
+## alpha       1.374   5.355 2.296e-04  0.9009  2.097
+## beta        1.234   3.398 3.949e-03  0.6341  2.401
 ## 
 ## Chi2 error levels in percent:
 ##          err.min n.optim df
@@ -638,10 +645,10 @@ plot(m.L2.DFOP)
 
summary(m.L2.DFOP, data = FALSE)
 
-
## mkin version:    0.9.35 
-## R version:       3.1.2 
-## Date of fit:     Sat Feb 21 14:44:57 2015 
-## Date of summary: Sat Feb 21 14:44:57 2015 
+
## mkin version:    0.9.36 
+## R version:       3.2.1 
+## Date of fit:     Sun Jun 21 01:48:02 2015 
+## Date of summary: Sun Jun 21 01:48:02 2015 
 ## 
 ## Equations:
 ## d_parent = - ((k1 * g * exp(-k1 * time) + k2 * (1 - g) * exp(-k2 *
@@ -650,7 +657,7 @@ plot(m.L2.DFOP)
 ## 
 ## Model predictions using solution type analytical 
 ## 
-## Fitted with method Port using 336 model solutions performed in 0.856 s
+## Fitted with method Port using 336 model solutions performed in 0.835 s
 ## 
 ## Weighting: none
 ## 
@@ -671,12 +678,12 @@ plot(m.L2.DFOP)
 ## Fixed parameter values:
 ## None
 ## 
-## Optimised, transformed parameters:
-##          Estimate Std. Error Lower Upper t value Pr(>|t|) Pr(>t)
-## parent_0  93.9500         NA    NA    NA      NA       NA     NA
-## log_k1     3.1210         NA    NA    NA      NA       NA     NA
-## log_k2    -1.0880         NA    NA    NA      NA       NA     NA
-## g_ilr     -0.2821         NA    NA    NA      NA       NA     NA
+## Optimised, transformed parameters with symmetric confidence intervals:
+##          Estimate Std. Error Lower Upper
+## parent_0  93.9500         NA    NA    NA
+## log_k1     3.1210         NA    NA    NA
+## log_k2    -1.0880         NA    NA    NA
+## g_ilr     -0.2821         NA    NA    NA
 ## 
 ## Parameter correlation:
 ## Could not estimate covariance matrix; singular system:
@@ -684,11 +691,14 @@ plot(m.L2.DFOP)
 ## Residual standard error: 1.732 on 8 degrees of freedom
 ## 
 ## Backtransformed parameters:
-##          Estimate Lower Upper
-## parent_0  93.9500    NA    NA
-## k1        22.6700    NA    NA
-## k2         0.3369    NA    NA
-## g          0.4016    NA    NA
+##   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  93.9500      NA     NA    NA    NA
+## k1        22.6700      NA     NA    NA    NA
+## k2         0.3369      NA     NA    NA    NA
+## g          0.4016      NA     NA    NA    NA
 ## 
 ## Chi2 error levels in percent:
 ##          err.min n.optim df
@@ -727,17 +737,17 @@ plot(m.L3.SFO)
 
summary(m.L3.SFO)
 
-
## mkin version:    0.9.35 
-## R version:       3.1.2 
-## Date of fit:     Sat Feb 21 14:44:57 2015 
-## Date of summary: Sat Feb 21 14:44:57 2015 
+
## mkin version:    0.9.36 
+## R version:       3.2.1 
+## Date of fit:     Sun Jun 21 01:48:03 2015 
+## Date of summary: Sun Jun 21 01:48:03 2015 
 ## 
 ## Equations:
 ## d_parent = - k_parent_sink * parent
 ## 
 ## Model predictions using solution type analytical 
 ## 
-## Fitted with method Port using 43 model solutions performed in 0.109 s
+## Fitted with method Port using 43 model solutions performed in 0.104 s
 ## 
 ## Weighting: none
 ## 
@@ -754,13 +764,10 @@ plot(m.L3.SFO)
 ## Fixed parameter values:
 ## None
 ## 
-## Optimised, transformed parameters:
-##                   Estimate Std. Error  Lower Upper t value  Pr(>|t|)
-## parent_0            74.870     8.4570 54.180 95.57   8.853 1.155e-04
-## log_k_parent_sink   -3.678     0.3261 -4.476 -2.88 -11.280 2.903e-05
-##                      Pr(>t)
-## parent_0          5.776e-05
-## log_k_parent_sink 1.451e-05
+## Optimised, transformed parameters with symmetric confidence intervals:
+##                   Estimate Std. Error  Lower Upper
+## parent_0            74.870     8.4570 54.180 95.57
+## log_k_parent_sink   -3.678     0.3261 -4.476 -2.88
 ## 
 ## Parameter correlation:
 ##                   parent_0 log_k_parent_sink
@@ -770,9 +777,12 @@ plot(m.L3.SFO)
 ## Residual standard error: 12.91 on 6 degrees of freedom
 ## 
 ## Backtransformed parameters:
-##               Estimate    Lower    Upper
-## parent_0      74.87000 54.18000 95.57000
-## k_parent_sink  0.02527  0.01138  0.05612
+##   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      74.87000   8.853 5.776e-05 54.18000 95.57000
+## k_parent_sink  0.02527   3.067 1.102e-02  0.01138  0.05612
 ## 
 ## Chi2 error levels in percent:
 ##          err.min n.optim df
@@ -813,17 +823,17 @@ plot(m.L3.FOMC)
 
summary(m.L3.FOMC, data = FALSE)
 
-
## mkin version:    0.9.35 
-## R version:       3.1.2 
-## Date of fit:     Sat Feb 21 14:44:58 2015 
-## Date of summary: Sat Feb 21 14:44:58 2015 
+
## mkin version:    0.9.36 
+## R version:       3.2.1 
+## Date of fit:     Sun Jun 21 01:48:03 2015 
+## Date of summary: Sun Jun 21 01:48:03 2015 
 ## 
 ## Equations:
-## d_parent = - (alpha/beta) * ((time/beta) + 1)^-1 * parent
+## d_parent = - (alpha/beta) * 1/((time/beta) + 1) * parent
 ## 
 ## Model predictions using solution type analytical 
 ## 
-## Fitted with method Port using 83 model solutions performed in 0.203 s
+## Fitted with method Port using 83 model solutions performed in 0.196 s
 ## 
 ## Weighting: none
 ## 
@@ -842,11 +852,11 @@ plot(m.L3.FOMC)
 ## Fixed parameter values:
 ## None
 ## 
-## Optimised, transformed parameters:
-##           Estimate Std. Error   Lower    Upper t value  Pr(>|t|)    Pr(>t)
-## parent_0   96.9700     4.5500 85.2800 108.7000  21.310 4.216e-06 2.108e-06
-## log_alpha  -0.8619     0.1704 -1.3000  -0.4238  -5.057 3.911e-03 1.955e-03
-## log_beta    0.6193     0.4744 -0.6003   1.8390   1.305 2.486e-01 1.243e-01
+## Optimised, transformed parameters with symmetric confidence intervals:
+##           Estimate Std. Error   Lower    Upper
+## parent_0   96.9700     4.5500 85.2800 108.7000
+## log_alpha  -0.8619     0.1704 -1.3000  -0.4238
+## log_beta    0.6193     0.4744 -0.6003   1.8390
 ## 
 ## Parameter correlation:
 ##           parent_0 log_alpha log_beta
@@ -857,10 +867,13 @@ plot(m.L3.FOMC)
 ## Residual standard error: 4.572 on 5 degrees of freedom
 ## 
 ## Backtransformed parameters:
-##          Estimate   Lower    Upper
-## parent_0  96.9700 85.2800 108.7000
-## alpha      0.4224  0.2725   0.6546
-## beta       1.8580  0.5487   6.2890
+##   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  96.9700  21.310 2.108e-06 85.2800 108.7000
+## alpha      0.4224   5.867 1.020e-03  0.2725   0.6546
+## beta       1.8580   2.108 4.444e-02  0.5487   6.2890
 ## 
 ## Chi2 error levels in percent:
 ##          err.min n.optim df
@@ -886,10 +899,10 @@ plot(m.L3.DFOP)
 
summary(m.L3.DFOP, data = FALSE)
 
-
## mkin version:    0.9.35 
-## R version:       3.1.2 
-## Date of fit:     Sat Feb 21 14:44:58 2015 
-## Date of summary: Sat Feb 21 14:44:58 2015 
+
## mkin version:    0.9.36 
+## R version:       3.2.1 
+## Date of fit:     Sun Jun 21 01:48:03 2015 
+## Date of summary: Sun Jun 21 01:48:03 2015 
 ## 
 ## Equations:
 ## d_parent = - ((k1 * g * exp(-k1 * time) + k2 * (1 - g) * exp(-k2 *
@@ -898,7 +911,7 @@ plot(m.L3.DFOP)
 ## 
 ## Model predictions using solution type analytical 
 ## 
-## Fitted with method Port using 137 model solutions performed in 0.346 s
+## Fitted with method Port using 137 model solutions performed in 0.35 s
 ## 
 ## Weighting: none
 ## 
@@ -919,12 +932,12 @@ plot(m.L3.DFOP)
 ## Fixed parameter values:
 ## None
 ## 
-## Optimised, transformed parameters:
-##          Estimate Std. Error   Lower     Upper t value  Pr(>|t|)    Pr(>t)
-## parent_0  97.7500    1.43800 93.7500 101.70000  67.970 2.808e-07 1.404e-07
-## log_k1    -0.6612    0.13340 -1.0310  -0.29100  -4.958 7.715e-03 3.858e-03
-## log_k2    -4.2860    0.05902 -4.4500  -4.12200 -72.620 2.155e-07 1.077e-07
-## g_ilr     -0.1229    0.05121 -0.2651   0.01925  -2.401 7.431e-02 3.716e-02
+## Optimised, transformed parameters with symmetric confidence intervals:
+##          Estimate Std. Error   Lower     Upper
+## parent_0  97.7500    1.43800 93.7500 101.70000
+## log_k1    -0.6612    0.13340 -1.0310  -0.29100
+## log_k2    -4.2860    0.05902 -4.4500  -4.12200
+## g_ilr     -0.1229    0.05121 -0.2651   0.01925
 ## 
 ## Parameter correlation:
 ##          parent_0  log_k1   log_k2   g_ilr
@@ -936,11 +949,14 @@ plot(m.L3.DFOP)
 ## Residual standard error: 1.439 on 4 degrees of freedom
 ## 
 ## Backtransformed parameters:
-##          Estimate    Lower     Upper
-## parent_0 97.75000 93.75000 101.70000
-## k1        0.51620  0.35650   0.74750
-## k2        0.01376  0.01168   0.01621
-## g         0.45660  0.40730   0.50680
+##   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 97.75000  67.970 1.404e-07 93.75000 101.70000
+## k1        0.51620   7.499 8.460e-04  0.35650   0.74750
+## k2        0.01376  16.940 3.557e-05  0.01168   0.01621
+## g         0.45660  25.410 7.121e-06  0.40730   0.50680
 ## 
 ## Chi2 error levels in percent:
 ##          err.min n.optim df
@@ -984,17 +1000,17 @@ plot(m.L4.SFO)
 
summary(m.L4.SFO, data = FALSE)
 
-
## mkin version:    0.9.35 
-## R version:       3.1.2 
-## Date of fit:     Sat Feb 21 14:44:58 2015 
-## Date of summary: Sat Feb 21 14:44:58 2015 
+
## mkin version:    0.9.36 
+## R version:       3.2.1 
+## Date of fit:     Sun Jun 21 01:48:04 2015 
+## Date of summary: Sun Jun 21 01:48:04 2015 
 ## 
 ## Equations:
 ## d_parent = - k_parent_sink * parent
 ## 
 ## Model predictions using solution type analytical 
 ## 
-## Fitted with method Port using 46 model solutions performed in 0.109 s
+## Fitted with method Port using 46 model solutions performed in 0.105 s
 ## 
 ## Weighting: none
 ## 
@@ -1011,13 +1027,10 @@ plot(m.L4.SFO)
 ## Fixed parameter values:
 ## None
 ## 
-## Optimised, transformed parameters:
-##                   Estimate Std. Error  Lower   Upper t value  Pr(>|t|)
-## parent_0             96.44    1.94900 91.670 101.200   49.49 4.566e-09
-## log_k_parent_sink    -5.03    0.07999 -5.225  -4.834  -62.88 1.088e-09
-##                      Pr(>t)
-## parent_0          2.283e-09
-## log_k_parent_sink 5.438e-10
+## Optimised, transformed parameters with symmetric confidence intervals:
+##                   Estimate Std. Error  Lower   Upper
+## parent_0             96.44    1.94900 91.670 101.200
+## log_k_parent_sink    -5.03    0.07999 -5.225  -4.834
 ## 
 ## Parameter correlation:
 ##                   parent_0 log_k_parent_sink
@@ -1027,9 +1040,12 @@ plot(m.L4.SFO)
 ## Residual standard error: 3.651 on 6 degrees of freedom
 ## 
 ## Backtransformed parameters:
-##                Estimate     Lower     Upper
-## parent_0      96.440000 91.670000 1.012e+02
-## k_parent_sink  0.006541  0.005378 7.955e-03
+##   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      96.440000   49.49 2.283e-09 91.670000 1.012e+02
+## k_parent_sink  0.006541   12.50 8.008e-06  0.005378 7.955e-03
 ## 
 ## Chi2 error levels in percent:
 ##          err.min n.optim df
@@ -1059,17 +1075,17 @@ plot(m.L4.FOMC)
 
summary(m.L4.FOMC, data = FALSE)
 
-
## mkin version:    0.9.35 
-## R version:       3.1.2 
-## Date of fit:     Sat Feb 21 14:44:58 2015 
-## Date of summary: Sat Feb 21 14:44:58 2015 
+
## mkin version:    0.9.36 
+## R version:       3.2.1 
+## Date of fit:     Sun Jun 21 01:48:04 2015 
+## Date of summary: Sun Jun 21 01:48:04 2015 
 ## 
 ## Equations:
-## d_parent = - (alpha/beta) * ((time/beta) + 1)^-1 * parent
+## d_parent = - (alpha/beta) * 1/((time/beta) + 1) * parent
 ## 
 ## Model predictions using solution type analytical 
 ## 
-## Fitted with method Port using 66 model solutions performed in 0.161 s
+## Fitted with method Port using 66 model solutions performed in 0.159 s
 ## 
 ## Weighting: none
 ## 
@@ -1088,11 +1104,11 @@ plot(m.L4.FOMC)
 ## Fixed parameter values:
 ## None
 ## 
-## Optimised, transformed parameters:
-##           Estimate Std. Error  Lower    Upper t value  Pr(>|t|)    Pr(>t)
-## parent_0   99.1400     1.6800 94.820 103.5000 59.0200 2.643e-08 1.322e-08
-## log_alpha  -0.3506     0.3725 -1.308   0.6068 -0.9414 3.897e-01 1.949e-01
-## log_beta    4.1740     0.5635  2.726   5.6230  7.4070 7.059e-04 3.530e-04
+## Optimised, transformed parameters with symmetric confidence intervals:
+##           Estimate Std. Error  Lower    Upper
+## parent_0   99.1400     1.6800 94.820 103.5000
+## log_alpha  -0.3506     0.3725 -1.308   0.6068
+## log_beta    4.1740     0.5635  2.726   5.6230
 ## 
 ## Parameter correlation:
 ##           parent_0 log_alpha log_beta
@@ -1103,10 +1119,13 @@ plot(m.L4.FOMC)
 ## Residual standard error: 2.315 on 5 degrees of freedom
 ## 
 ## Backtransformed parameters:
-##          Estimate   Lower   Upper
-## parent_0  99.1400 94.8200 103.500
-## alpha      0.7042  0.2703   1.835
-## beta      64.9800 15.2600 276.600
+##   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  99.1400  59.020 1.322e-08 94.8200 103.500
+## alpha      0.7042   2.685 2.178e-02  0.2703   1.835
+## beta      64.9800   1.775 6.807e-02 15.2600 276.600
 ## 
 ## Chi2 error levels in percent:
 ##          err.min n.optim df
diff --git a/vignettes/FOCUS_Z.Rnw b/vignettes/FOCUS_Z.Rnw
index 5e2e0251..1df0ee9c 100644
--- a/vignettes/FOCUS_Z.Rnw
+++ b/vignettes/FOCUS_Z.Rnw
@@ -265,11 +265,13 @@ summary(m.Z.mkin.5a, data = FALSE)$bpar
 @
 
 A graphical representation of the confidence intervals can finally be obtained.
+
 <>=
 mkinparplot(m.Z.mkin.5a)
 @
 
 The endpoints obtained with this model are
+
 <>=
 endpoints(m.Z.mkin.5a)
 @
diff --git a/vignettes/FOCUS_Z.pdf b/vignettes/FOCUS_Z.pdf
index 3174a23a..36f3dc14 100644
Binary files a/vignettes/FOCUS_Z.pdf and b/vignettes/FOCUS_Z.pdf differ
diff --git a/vignettes/compiled_models.html b/vignettes/compiled_models.html
index 2f2a6edb..e6f21b09 100644
--- a/vignettes/compiled_models.html
+++ b/vignettes/compiled_models.html
@@ -77,7 +77,7 @@ img {
 -->
 

Benchmark for a model that can also be solved with Eigenvalues

-

This evaluation is taken from the example section of mkinfit. When using an mkin version greater than 0.9-36 and the ccSolve package is installed and functional, you will get a message that the model is being compiled when defining a model using mkinmod.

+

This evaluation is taken from the example section of mkinfit. When using an mkin version equal to or greater than 0.9-36 and a compiler (gcc) is installed, you will see a message that the model is being compiled from autogenerated C code when defining a model using mkinmod.

library("mkin")
 SFO_SFO <- mkinmod(
   parent = list(type = "SFO", to = "m1", sink = TRUE),
@@ -94,20 +94,20 @@ mb.1 <- microbenchmark(
 smb.1 <- summary(mb.1)[-1]
 rownames(smb.1) <- c("deSolve, not compiled", "Eigenvalue based", "deSolve, compiled")
 print(smb.1)
-
##                             min        lq      mean    median        uq
-## deSolve, not compiled 6192.0125 6195.3470 6211.0309 6198.6816 6220.5401
-## Eigenvalue based       956.7604 1008.7224 1026.2572 1060.6844 1061.0055
-## deSolve, compiled      869.6880  871.9315  883.4929  874.1751  890.3953
+
##                            min        lq      mean    median        uq
+## deSolve, not compiled 4969.585 5033.7311 5092.7389 5097.8773 5154.3160
+## Eigenvalue based       868.731  891.7239  909.6449  914.7169  930.1018
+## deSolve, compiled     4935.049 4935.4796 4968.2150 4935.9097 4984.7978
 ##                             max neval
-## deSolve, not compiled 6242.3986     3
-## Eigenvalue based      1061.3266     3
-## deSolve, compiled      906.6155     3
+## deSolve, not compiled 5210.7547 3 +## Eigenvalue based 945.4867 3 +## deSolve, compiled 5033.6858 3

We see that using the compiled model is almost a factor of 8 faster than using the R version with the default ode solver, and it is even faster than the Eigenvalue based solution implemented in R which does not need iterative solution of the ODEs:

smb.1["median"]/smb.1["deSolve, compiled", "median"]
-
##                         median
-## deSolve, not compiled 7.120877
-## Eigenvalue based      1.205328
-## deSolve, compiled     1.000000
+
##                          median
+## deSolve, not compiled 1.0328141
+## Eigenvalue based      0.1853188
+## deSolve, compiled     1.0000000

Benchmark for a model that can not be solved with Eigenvalues

@@ -124,16 +124,16 @@ smb.2 <- summary(mb.2)[-1] rownames(smb.2) <- c("deSolve, not compiled", "deSolve, compiled") print(smb.2)
##                             min        lq      mean    median        uq
-## deSolve, not compiled 13.297283 13.427702 13.481155 13.558121 13.573092
-## deSolve, compiled      1.486926  1.526887  1.546851  1.566848  1.576813
-##                             max neval
-## deSolve, not compiled 13.588063     3
-## deSolve, compiled      1.586778     3
+## deSolve, not compiled 11.745276 11.754288 11.820726 11.763300 11.858451 +## deSolve, compiled 1.385829 1.386407 1.400841 1.386985 1.408347 +## max neval +## deSolve, not compiled 11.95360 3 +## deSolve, compiled 1.42971 3
smb.2["median"]/smb.2["deSolve, compiled", "median"]
-
##                         median
-## deSolve, not compiled 8.653119
-## deSolve, compiled     1.000000
-

Here we get a performance benefit of more than a factor of 8 using the version of the differential equation model compiled from C code using the ccSolve package!

+
##                       median
+## deSolve, not compiled 8.4812
+## deSolve, compiled     1.0000
+

Here we get a performance benefit of more than a factor of 10 using the version of the differential equation model compiled from C code using the inline package!

diff --git a/vignettes/mkin.pdf b/vignettes/mkin.pdf index 8f1e0884..5786c5bf 100644 Binary files a/vignettes/mkin.pdf and b/vignettes/mkin.pdf differ -- cgit v1.2.1