- Testing of the new fitting process with transformed parameters shows that it is

  less stable than the way of fitting used in mkin 0.8. Presumably this is due
  to the frequent presence of products of two parameters (formation fractions and kinetic
  rate constants) in the differential equations.
- mkinmod documentation was adapted to the new version
- Introduce vim fold markers to the mkinmod source for convenience




git-svn-id: svn+ssh://svn.r-forge.r-project.org/svnroot/kinfit/pkg/mkin@24 edb9625f-4e0d-4859-8d74-9fd3b1da38cb

1 files changed, 37 insertions, 40 deletions

diff --git a/R/mkinmod.R b/R/mkinmod.R
index d4ff75b..37f96ce 100644
--- a/R/mkinmod.R
+++ b/R/mkinmod.R
@@ -1,6 +1,6 @@
 # $Id$

 

-# Copyright (C) 2010-2012 Johannes Ranke

+# Copyright (C) 2010-2012 Johannes Ranke#{{{

 # Contact: jranke@uni-bremen.de

 

 # This file is part of the R package mkin

@@ -16,20 +16,23 @@
 # details.

 

 # You should have received a copy of the GNU General Public License along with

-# this program. If not, see <http://www.gnu.org/licenses/>

+# this program. If not, see <http://www.gnu.org/licenses/>#}}}

 

 mkinmod <- function(...)

 {

   spec <- list(...)

   obs_vars <- names(spec)

 

-  # The returned model will be a list of character vectors, containing

+  # The returned model will be a list of character vectors, containing#{{{

   # differential equations, parameter names and a mapping from model variables

   # to observed variables. If possible, a matrix representation of the 

   # differential equations is included

   parms <- vector()

   diffs <- vector()

-  map <- list()

+  map <- list()#}}}

+

+  # Give a warning when a model with time dependent degradation uses formation#{{{

+  # fractions

   if(spec[[1]]$type %in% c("FOMC", "DFOP", "HS")) {

     mat = FALSE 

     if(!is.null(spec[[1]]$to)) {

@@ -40,12 +43,13 @@ mkinmod <- function(...)
         sep="\n")

       warning(message)

     } else message <- "ok"

-  } else mat = TRUE

+  } else mat = TRUE#}}}

 

-  # Establish list of differential equations

-  # as well as map from observed compartments to differential equations

+  # Establish list of differential equations as well as map from observed#{{{

+  # compartments to differential equations

   for (varname in obs_vars)

   {

+    # Check the type component of the compartment specification#{{{

     if(is.null(spec[[varname]]$type)) stop(

       "Every argument to mkinmod must be a list containing a type component")

     if(!spec[[varname]]$type %in% c("SFO", "FOMC", "DFOP", "HS", "SFORB")) stop(

@@ -53,9 +57,8 @@ mkinmod <- function(...)
     if(spec[[varname]]$type %in% c("FOMC", "DFOP", "HS") & match(varname, obs_vars) != 1) {

         stop(paste("Types FOMC, DFOP and HS are only implemented for the first compartment,", 

                    "which is assumed to be the source compartment"))

-    }

-

-    # New (sub)compartments (boxes) needed for the model type

+    }#}}}

+    # New (sub)compartments (boxes) needed for the model type#{{{

     new_boxes <- switch(spec[[varname]]$type,

       SFO = varname,

       FOMC = varname,

@@ -64,51 +67,44 @@ mkinmod <- function(...)
       SFORB = paste(varname, c("free", "bound"), sep="_")

     )

     map[[varname]] <- new_boxes

-    names(map[[varname]]) <- rep(spec[[varname]]$type, length(new_boxes))

-

-    # Start a new differential equation for each new box

+    names(map[[varname]]) <- rep(spec[[varname]]$type, length(new_boxes))#}}}

+    # Start a new differential equation for each new box#{{{

     new_diffs <- paste("d_", new_boxes, " =", sep="")

     names(new_diffs) <- new_boxes

-    diffs <- c(diffs, new_diffs)

-  }

+    diffs <- c(diffs, new_diffs)#}}}

+  }#}}}

 

-  # Create content of differential equations and build parameter list

+  # Create content of differential equations and build parameter list#{{{

   for (varname in obs_vars)

   {

-    # Add first-order term to first (or only) box for SFO and SFORB

+    # Add first-order term to first (or only) box for SFO and SFORB#{{{

     box_1 = map[[varname]][[1]] # This is the only box unless type is SFORB

     if(spec[[varname]]$type %in% c("SFO", "SFORB")) {

       k_compound <- paste("k", box_1, sep="_")

       origin_term <- paste(k_compound, "*", box_1)

       parms <- c(parms, k_compound)

-    }

-

-

-    # Construct and add FOMC term and add FOMC parameters if needed

+    }#}}}

+    # Construct and add FOMC term and add FOMC parameters if needed#{{{

     if(spec[[varname]]$type == "FOMC") {

       # From p. 53 of the FOCUS kinetics report

       origin_term <- paste("(alpha/beta) * ((time/beta) + 1)^-1 *", box_1)

       parms <- c(parms, "alpha", "beta")

-    } 

-

-    # Construct and add DFOP term and add DFOP parameters if needed

+    } #}}}

+    # Construct and add DFOP term and add DFOP parameters if needed#{{{

     if(spec[[varname]]$type == "DFOP") {

       # From p. 57 of the FOCUS kinetics report

       origin_term <- paste("((k1 * g * exp(-k1 * time) + k2 * (1 - g) * exp(-k2 * time)) / (g * exp(-k1 * time) + (1 - g) * exp(-k2 * time))) *", box_1)

       parms <- c(parms, "k1", "k2", "g")

-    } 

-

-    # Construct and add HS term and add HS parameters if needed

+    } #}}}

+    # Construct and add HS term and add HS parameters if needed#{{{

     if(spec[[varname]]$type == "HS") {

       # From p. 55 of the FOCUS kinetics report

       origin_term <- paste("ifelse(time <= tb, k1, k2)", "*", box_1)

       parms <- c(parms, "k1", "k2", "tb")

-    } 

-

-    # Add origin decline term to box 1 (usually the only box, unless type is SFORB)

-    diffs[[box_1]] <- paste(diffs[[box_1]], "-", origin_term)

-   

-    # Add reversible binding for SFORB models

+    } #}}}

+    # Add origin decline term to box 1 (usually the only box, unless type is SFORB)#{{{

+    diffs[[box_1]] <- paste(diffs[[box_1]], "-", origin_term)#}}}

+    # Add reversible binding for SFORB models#{{{

     if(spec[[varname]]$type == "SFORB") {

       box_2 = map[[varname]][[2]]

       k_free_bound <- paste("k", varname, "free", "bound", sep="_")

@@ -120,9 +116,8 @@ mkinmod <- function(...)
         k_bound_free, "*", box_2)

       diffs[[box_2]] <- paste(diffs[[box_2]], reversible_binding_term_2)

       parms <- c(parms, k_free_bound, k_bound_free)

-    } 

-

-    # Transfer between compartments

+    } #}}}

+    # Transfer between compartments#{{{

     to <- spec[[varname]]$to

     if(!is.null(to)) {

       # Name of box from which transfer takes place

@@ -138,11 +133,12 @@ mkinmod <- function(...)
             fraction_to_target, "*", origin_term)

         parms <- c(parms, fraction_to_target)

       }

-    }

-  }

+    }#}}}

+  }#}}}

+

   model <- list(diffs = diffs, parms = parms, map = map)

 

-  # Create coefficient matrix if appropriate

+  # Create coefficient matrix if appropriate#{{{

   if (mat) {

     boxes <- names(diffs)

     n <- length(boxes)

@@ -165,8 +161,9 @@ mkinmod <- function(...)
       }

     }

     model$coefmat <- m

-  }

+  }#}}}

 

   class(model) <- "mkinmod"

   return(model)

 }

+# vim: set foldmethod=marker: