From c1144753adfa0809003085009ebd85f8af9beda8 Mon Sep 17 00:00:00 2001
From: jranke <jranke@edb9625f-4e0d-4859-8d74-9fd3b1da38cb>
Date: Tue, 10 Apr 2012 21:50:22 +0000
Subject: - Fitting and summaries now work with the new parameter
 transformations. - The SFORB models with metabolites is broken (see TODO) -
 Moved the vignette to the location recommended since R 2.14 - Added the
 missing documentation - Commented out the schaefer_complex_case test, as this
 version of   mkin is not able to fit a model without sink and therefore  
 mkin estimated parameters are quite different

git-svn-id: svn+ssh://svn.r-forge.r-project.org/svnroot/kinfit/pkg/mkin@22 edb9625f-4e0d-4859-8d74-9fd3b1da38cb
---
 R/transform_odeparms.R | 70 ++++++++++++++++++++++++++++++++++----------------
 1 file changed, 48 insertions(+), 22 deletions(-)

(limited to 'R/transform_odeparms.R')

diff --git a/R/transform_odeparms.R b/R/transform_odeparms.R
index 4b8bfd14..8b39874e 100644
--- a/R/transform_odeparms.R
+++ b/R/transform_odeparms.R
@@ -1,46 +1,72 @@
-transform_odeparms <- function(odeparms, mod_vars) {
-  # Transform rate constants and formation fractions
-  transparms <- odeparms
-  # Exponential transformation for rate constants
-  index_k <- grep("^k_", names(odeparms))
+transform_odeparms <- function(parms, mod_vars) {
+  # Set up container for transformed parameters
+  transparms <- parms
+
+  # Log transformation for rate constants
+  index_k <- grep("^k_", names(transparms))
   if (length(index_k) > 0) {
-    transparms[index_k] <- exp(odeparms[index_k])
+    transparms[index_k] <- log(parms[index_k])
   }
 
-  # Go through state variables and apply inverse isotropic logratio transformation
+  # Go through state variables and apply isotropic logratio transformation
   for (box in mod_vars) {
-    indices_f <- grep(paste("^f", box, sep = "_"), names(odeparms))
-    f_names <- grep(paste("^f", box, sep = "_"), names(odeparms), value = TRUE)
+    indices_f <- grep(paste("^f", box, sep = "_"), names(parms))
+    f_names <- grep(paste("^f", box, sep = "_"), names(parms), value = TRUE)
     n_paths <- length(indices_f)
     if (n_paths > 0) {
-      f <- invilr(odeparms[indices_f])[1:n_paths] # We do not need the last component
-      names(f) <- f_names
-      transparms[indices_f] <- f
+      f <- parms[indices_f]
+      trans_f <- ilr(c(f, 1 - sum(f)))
+      names(trans_f) <- f_names
+      transparms[indices_f] <- trans_f
     }
   }
+
+  # Transform parameters also for FOMC, DFOP and HS models
+  for (pname in c("alpha", "beta", "k1", "k2", "tb")) {
+    if (!is.na(parms[pname])) {
+      transparms[pname] <- log(parms[pname])
+    } 
+  }
+  if (!is.na(parms["g"])) {
+    g <- parms["g"]
+    transparms["g"] <- ilr(c(g, 1- g))
+  }
+
   return(transparms)
 }
 
 backtransform_odeparms <- function(transparms, mod_vars) {
-  # Transform rate constants and formation fractions
-  odeparms <- transparms
-  # Log transformation for rate constants
-  index_k <- grep("^k_", names(transparms))
+  # Set up container for backtransformed parameters
+  parms <- transparms
+
+  # Exponential transformation for rate constants
+  index_k <- grep("^k_", names(parms))
   if (length(index_k) > 0) {
-    odeparms[index_k] <- log(transparms[index_k])
+    parms[index_k] <- exp(transparms[index_k])
   }
 
-  # Go through state variables and apply isotropic logratio transformation
+  # Go through state variables and apply inverse isotropic logratio transformation
   for (box in mod_vars) {
     indices_f <- grep(paste("^f", box, sep = "_"), names(transparms))
     f_names <- grep(paste("^f", box, sep = "_"), names(transparms), value = TRUE)
     n_paths <- length(indices_f)
     if (n_paths > 0) {
-      trans_f <- transparms[indices_f]
-      f <- ilr(c(trans_f, 1 - sum(trans_f)))
+      f <- invilr(transparms[indices_f])[1:n_paths] # We do not need the last component
       names(f) <- f_names
-      odeparms[indices_f] <- f
+      parms[indices_f] <- f
     }
   }
-  return(odeparms)
+
+  # Transform parameters also for DFOP and HS models
+  for (pname in c("alpha", "beta", "k1", "k2", "tb")) {
+    if (!is.na(transparms[pname])) {
+      parms[pname] <- exp(transparms[pname])
+    } 
+  }
+  if (!is.na(transparms["g"])) {
+    g <- transparms["g"]
+    parms["g"] <- invilr(g)[1]
+  }
+
+  return(parms)
 }
-- 
cgit v1.2.1