aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--R/mkinfit.R3
-rw-r--r--R/mkinpredict.R3
-rw-r--r--inst/unitTests/runit.mkinerrmin.R12
-rw-r--r--man/mkinfit.Rd21
-rw-r--r--man/mkinpredict.Rd24
5 files changed, 56 insertions, 7 deletions
diff --git a/R/mkinfit.R b/R/mkinfit.R
index b0df6e8c..38f76674 100644
--- a/R/mkinfit.R
+++ b/R/mkinfit.R
@@ -28,6 +28,7 @@ mkinfit <- function(mkinmod, observed,
fixed_initials = names(mkinmod$diffs)[-1],
solution_type = "auto",
method.ode = "lsoda",
+ use_compiled = "auto",
method.modFit = c("Port", "Marq", "SANN", "Nelder-Mead", "BFGS", "CG", "L-BFGS-B"),
maxit.modFit = "auto",
control.modFit = list(),
@@ -263,6 +264,7 @@ mkinfit <- function(mkinmod, observed,
out <- mkinpredict(mkinmod, parms,
odeini, outtimes,
solution_type = solution_type,
+ use_compiled = use_compiled,
method.ode = method.ode,
atol = atol, rtol = rtol, ...)
@@ -282,6 +284,7 @@ mkinfit <- function(mkinmod, observed,
out_plot <- mkinpredict(mkinmod, parms,
odeini, outtimes_plot,
solution_type = solution_type,
+ use_compiled = use_compiled,
method.ode = method.ode,
atol = atol, rtol = rtol, ...)
diff --git a/R/mkinpredict.R b/R/mkinpredict.R
index bed82441..864eefbe 100644
--- a/R/mkinpredict.R
+++ b/R/mkinpredict.R
@@ -19,6 +19,7 @@
mkinpredict <- function(mkinmod, odeparms, odeini,
outtimes, solution_type = "deSolve",
+ use_compiled = "auto",
method.ode = "lsoda", atol = 1e-8, rtol = 1e-10,
map_output = TRUE, ...) {
@@ -87,7 +88,7 @@ mkinpredict <- function(mkinmod, odeparms, odeini,
names(out) <- c("time", mod_vars)
}
if (solution_type == "deSolve") {
- if (!is.null(mkinmod$compiled)) {
+ if (!is.null(mkinmod$compiled) & use_compiled[1] != FALSE) {
mkindiff <- mkinmod$compiled
} else {
mkindiff <- function(t, state, parms) {
diff --git a/inst/unitTests/runit.mkinerrmin.R b/inst/unitTests/runit.mkinerrmin.R
index 56a33ff9..96793a36 100644
--- a/inst/unitTests/runit.mkinerrmin.R
+++ b/inst/unitTests/runit.mkinerrmin.R
@@ -8,26 +8,34 @@ test.FOCUS_2006_D_SFO_SFO <- function()
m1 = list(type = "SFO"), use_of_ff = "max")
fit.1.e <- mkinfit(SFO_SFO.1, FOCUS_2006_D)
- fit.1.d <- mkinfit(SFO_SFO.1, solution_type = "deSolve", FOCUS_2006_D)
+ fit.1.d <- mkinfit(SFO_SFO.1, solution_type = "deSolve", use_compiled = FALSE, FOCUS_2006_D)
+ fit.1.dc <- mkinfit(SFO_SFO.1, solution_type = "deSolve", use_compiled = TRUE, FOCUS_2006_D)
fit.2.e <- mkinfit(SFO_SFO.2, FOCUS_2006_D)
- fit.2.d <- mkinfit(SFO_SFO.2, solution_type = "deSolve", FOCUS_2006_D)
+ fit.2.d <- mkinfit(SFO_SFO.2, solution_type = "deSolve", use_compiled = FALSE, FOCUS_2006_D)
+ fit.2.dc <- mkinfit(SFO_SFO.2, solution_type = "deSolve", use_compiled = TRUE, FOCUS_2006_D)
FOCUS_2006_D_results_schaefer07_means <- c(
parent_0 = 99.65, DT50_parent = 7.04, DT50_m1 = 131.34)
r.1.e <- c(fit.1.e$bparms.optim[[1]], endpoints(fit.1.e)$distimes[[1]])
r.1.d <- c(fit.1.d$bparms.optim[[1]], endpoints(fit.1.d)$distimes[[1]])
+ r.1.dc <- c(fit.1.dc$bparms.optim[[1]], endpoints(fit.1.dc)$distimes[[1]])
r.2.e <- c(fit.2.e$bparms.optim[[1]], endpoints(fit.2.e)$distimes[[1]])
r.2.d <- c(fit.2.d$bparms.optim[[1]], endpoints(fit.2.d)$distimes[[1]])
+ r.2.dc <- c(fit.2.dc$bparms.optim[[1]], endpoints(fit.2.dc)$distimes[[1]])
dev.1.e <- 100 * (r.1.e - FOCUS_2006_D_results_schaefer07_means)/r.1.e
checkIdentical(as.numeric(abs(dev.1.e)) < 1, rep(TRUE, 3))
dev.1.d <- 100 * (r.1.d - FOCUS_2006_D_results_schaefer07_means)/r.1.d
checkIdentical(as.numeric(abs(dev.1.d)) < 1, rep(TRUE, 3))
+ dev.1.dc <- 100 * (r.1.dc - FOCUS_2006_D_results_schaefer07_means)/r.1.dc
+ checkIdentical(as.numeric(abs(dev.1.dc)) < 1, rep(TRUE, 3))
dev.2.e <- 100 * (r.2.e - FOCUS_2006_D_results_schaefer07_means)/r.2.e
checkIdentical(as.numeric(abs(dev.2.e)) < 1, rep(TRUE, 3))
dev.2.d <- 100 * (r.2.d - FOCUS_2006_D_results_schaefer07_means)/r.2.d
checkIdentical(as.numeric(abs(dev.2.d)) < 1, rep(TRUE, 3))
+ dev.2.dc <- 100 * (r.2.dc - FOCUS_2006_D_results_schaefer07_means)/r.2.dc
+ checkIdentical(as.numeric(abs(dev.2.d)) < 1, rep(TRUE, 3))
round(mkinerrmin(fit.2.e), 4)
round(mkinerrmin(fit.2.d), 4)
diff --git a/man/mkinfit.Rd b/man/mkinfit.Rd
index 4b0ad1a7..4250bfda 100644
--- a/man/mkinfit.Rd
+++ b/man/mkinfit.Rd
@@ -23,6 +23,7 @@ mkinfit(mkinmod, observed,
fixed_parms = NULL, fixed_initials = names(mkinmod$diffs)[-1],
solution_type = "auto",
method.ode = "lsoda",
+ use_compiled = "auto",
method.modFit = c("Port", "Marq", "SANN", "Nelder-Mead", "BFGS", "CG", "L-BFGS-B"),
maxit.modFit = "auto",
control.modFit = list(),
@@ -102,6 +103,10 @@ mkinfit(mkinmod, observed,
\code{\link{ode}} in case the solution type is "deSolve". The default
"lsoda" is performant, but sometimes fails to converge.
}
+ \item{use_compiled}{
+ If set to \code{FALSE}, no compiled version of the \code{\link{mkinmod}}
+ model is used, even if is present.
+ }
\item{method.modFit}{
The optimisation method passed to \code{\link{modFit}}.
@@ -235,9 +240,18 @@ SFO_SFO <- mkinmod(
parent = list(type = "SFO", to = "m1", sink = TRUE),
m1 = list(type = "SFO"))
# Fit the model to the FOCUS example dataset D using defaults
-fit <- mkinfit(SFO_SFO, FOCUS_2006_D, solution_type = "eigen")
-fit.deSolve <- mkinfit(SFO_SFO, FOCUS_2006_D, solution_type = "deSolve")
-summary(fit)
+system.time(fit <- mkinfit(SFO_SFO, FOCUS_2006_D, solution_type = "eigen"))
+system.time(fit.deSolve <- mkinfit(SFO_SFO, FOCUS_2006_D, solution_type = "deSolve", use_compiled = FALSE))
+system.time(fit.deSolve.compiled <- mkinfit(SFO_SFO, FOCUS_2006_D, solution_type = "deSolve", use_compiled = TRUE))
+# The performance benefit of using the compiled version is immense compared with deSolve without compilation, and
+# about a factor of two compared with the eigenvalue based solution
+coef(fit)
+coef(fit.deSolve)
+coef(fit.deSolve.compiled)
+endpoints(fit)
+endpoints(fit.deSolve)
+# Check compiled model versions against other solutions
+
# Use stepwise fitting, using optimised parameters from parent only fit, FOMC
\dontrun{
@@ -289,6 +303,7 @@ f.w.man.irls <- mkinfit(SFO_SFO.ff, dw, err = "err.man",
reweight.method = "obs")
summary(f.w.man.irls)
}
+
}
\keyword{ models }
\keyword{ optimize }
diff --git a/man/mkinpredict.Rd b/man/mkinpredict.Rd
index 7d8979e4..7450a7ce 100644
--- a/man/mkinpredict.Rd
+++ b/man/mkinpredict.Rd
@@ -10,7 +10,8 @@
}
\usage{
mkinpredict(mkinmod, odeparms, odeini, outtimes, solution_type = "deSolve",
- method.ode = "lsoda", atol = 1e-08, rtol = 1e-10, map_output = TRUE, ...)
+ method.ode = "lsoda", use_compiled = "auto", atol = 1e-08, rtol = 1e-10,
+ map_output = TRUE, ...)
}
\arguments{
\item{mkinmod}{
@@ -41,6 +42,10 @@
\code{\link{ode}} in case the solution type is "deSolve". The default
"lsoda" is performant, but sometimes fails to converge.
}
+ \item{use_compiled}{
+ If set to \code{FALSE}, no compiled version of the \code{\link{mkinmod}}
+ model is used, even if is present.
+ }
\item{atol}{
Absolute error tolerance, passed to \code{\link{ode}}. Default is 1e-8,
lower than in \code{\link{lsoda}}.
@@ -71,6 +76,8 @@
mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20,
solution_type = "deSolve")
mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20,
+ solution_type = "deSolve", use_compiled = FALSE)
+ mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100), 0:20,
solution_type = "eigen")
@@ -92,5 +99,20 @@
mkinpredict(SFO, c(k_degradinol_sink = 0.3), c(degradinol = 100),
seq(0, 20, by = 0.01))[2001,]
+ # 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"))
+ system.time(
+ print(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),
+ solution_type = "eigen")[201,]))
+ system.time(
+ print(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),
+ solution_type = "deSolve")[201,]))
+ system.time(
+ print(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),
+ solution_type = "deSolve", use_compiled = FALSE)[201,]))
}
\keyword{ manip }

Contact - Imprint