From 5b28712ad11470293e7ab6746718ab3edcfe966b Mon Sep 17 00:00:00 2001
From: jranke <jranke@edb9625f-4e0d-4859-8d74-9fd3b1da38cb>
Date: Tue, 3 Jul 2012 15:46:14 +0000
Subject: - Add rtol argument in order to increase the precision of deSolve
 based   predictions for complex models - Improve default values for atol and
 set rtol default to 1e-10

git-svn-id: svn+ssh://svn.r-forge.r-project.org/svnroot/kinfit/pkg/mkin@49 edb9625f-4e0d-4859-8d74-9fd3b1da38cb
---
 R/mkinfit.R     | 7 ++++---
 R/mkinplot.R    | 2 +-
 R/mkinpredict.R | 3 ++-
 3 files changed, 7 insertions(+), 5 deletions(-)

(limited to 'R')

diff --git a/R/mkinfit.R b/R/mkinfit.R
index 3100b30..0215596 100644
--- a/R/mkinfit.R
+++ b/R/mkinfit.R
@@ -29,7 +29,7 @@ mkinfit <- function(mkinmod, observed,
   solution_type = "auto",
   plot = FALSE, quiet = FALSE,
   err = NULL, weight = "none", scaleVar = FALSE,
-  atol = 1e-6, n.outtimes = 100,
+  atol = 1e-8, rtol = 1e-10, n.outtimes = 100,
   trace_parms = FALSE,
   ...)
 {
@@ -130,7 +130,7 @@ mkinfit <- function(mkinmod, observed,
     parms <- backtransform_odeparms(odeparms, mod_vars)
 
     # Solve the system with current transformed parameter values
-    out <- mkinpredict(mkinmod, parms, odeini, outtimes, solution_type = solution_type, ...)
+    out <- mkinpredict(mkinmod, parms, odeini, outtimes, solution_type = solution_type, atol = atol, rtol = rtol, ...)
 
     assign("out_predicted", out, inherits=TRUE)
 
@@ -146,7 +146,7 @@ mkinfit <- function(mkinmod, observed,
         outtimes_plot = seq(min(observed$time), max(observed$time), length.out=100)
 
         out_plot <- mkinpredict(mkinmod, parms, odeini, outtimes_plot, 
-          solution_type = solution_type, ...)
+          solution_type = solution_type, atol = atol, rtol = rtol, ...)
 
         plot(0, type="n", 
           xlim = range(observed$time), ylim = range(observed$value, na.rm=TRUE),
@@ -198,6 +198,7 @@ mkinfit <- function(mkinmod, observed,
   data$variable <- ordered(data$variable, levels = obs_vars)
   fit$data <- data[order(data$variable, data$time), ]
   fit$atol <- atol
+  fit$rtol <- rtol
   fit$parms.all <- parms.all # Return all backtransformed parameters for summary
   fit$odeparms.final <- parms.all[mkinmod$parms] # Return ode parameters for further fitting
   fit$date <- date()
diff --git a/R/mkinplot.R b/R/mkinplot.R
index 7641d26..546c506 100644
--- a/R/mkinplot.R
+++ b/R/mkinplot.R
@@ -18,7 +18,7 @@ mkinplot <- function(fit, xlab = "Time", ylab = "Observed", xlim = range(fit$dat
   odeparms <- parms.all[odenames]
 
   out <- mkinpredict(fit$mkinmod, odeparms, odeini, outtimes, 
-          solution_type = solution_type, ...)
+          solution_type = solution_type, atol = fit$atol, rtol = fit$rtol, ...)
 
   # Plot the data and model output
   plot(0, type="n", 
diff --git a/R/mkinpredict.R b/R/mkinpredict.R
index 493c578..5170dc6 100644
--- a/R/mkinpredict.R
+++ b/R/mkinpredict.R
@@ -18,7 +18,7 @@
 # You should have received a copy of the GNU General Public License along with
 # this program. If not, see <http://www.gnu.org/licenses/>
 
-mkinpredict <- function(mkinmod, odeparms, odeini, outtimes, solution_type = "deSolve", map_output = TRUE, atol = 1e-6, ...) {
+mkinpredict <- function(mkinmod, odeparms, odeini, outtimes, solution_type = "deSolve", map_output = TRUE, atol = 1e-8, rtol = 1e-10, ...) {
 
   # Get the names of the state variables in the model
   mod_vars <- names(mkinmod$diffs)
@@ -92,6 +92,7 @@ mkinpredict <- function(mkinmod, odeparms, odeini, outtimes, solution_type = "de
       func = mkindiff, 
       parms = odeparms,
       atol = atol,
+      rtol = rtol,
       ...
     )
   }
-- 
cgit v1.2.1