aboutsummaryrefslogtreecommitdiff
path: root/R/mkinmod.R
diff options
context:
space:
mode:
authorjranke <jranke@edb9625f-4e0d-4859-8d74-9fd3b1da38cb>2010-05-18 12:58:38 +0000
committerjranke <jranke@edb9625f-4e0d-4859-8d74-9fd3b1da38cb>2010-05-18 12:58:38 +0000
commit16f5b1d3c0136413e92b2be0f20d365e92e9cd1c (patch)
treeedab7b27ba4253b9dd45520e504c54851f65f26f /R/mkinmod.R
parent30cbb4092f6d2d3beff5800603374a0d009ad770 (diff)
Much more complete version that was just submitted to CRAN.
git-svn-id: svn+ssh://svn.r-forge.r-project.org/svnroot/kinfit/pkg/mkin@9 edb9625f-4e0d-4859-8d74-9fd3b1da38cb
Diffstat (limited to 'R/mkinmod.R')
-rw-r--r--R/mkinmod.R66
1 files changed, 49 insertions, 17 deletions
diff --git a/R/mkinmod.R b/R/mkinmod.R
index 7e64f879..643e855e 100644
--- a/R/mkinmod.R
+++ b/R/mkinmod.R
@@ -1,38 +1,53 @@
-mkinmod <- function(spec = list(parent = list(type = "SFO", to = NA, sink = TRUE)))
+mkinmod <- function(...)
{
- if(!is.list(spec)) stop("spec must be a list of model specifications for each observed variable")
+ spec <- list(...)
+ obs_vars <- names(spec)
- # The returned model will be a list of two character vectors, containing
- # differential equations and parameter names
+ # The returned model will be a list of character vectors, containing
+ # differential equations, parameter names and a mapping from model variables
+ # to observed variables
parms <- vector()
diffs <- vector()
map <- list()
# Establish list of differential equations
- for (varname in names(spec))
+ for (varname in obs_vars)
{
- if(!spec[[varname]]$type %in% c("SFO", "SFORB")) stop("Available types are SFO and SFORB only")
+ 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", "SFORB")) stop("Available types are SFO, FOMC and SFORB only")
new_parms <- vector()
# New (sub)compartments (boxes) needed for the model type
new_boxes <- switch(spec[[varname]]$type,
SFO = varname,
+ FOMC = varname,
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
new_diffs <- paste("d_", new_boxes, " =", sep="")
# Construct terms for transfer to sink and add if appropriate
+ if(is.null(spec[[varname]]$sink)) spec[[varname]]$sink <- TRUE
if(spec[[varname]]$sink) {
- k_compound_sink <- paste("k", new_boxes[[1]], "sink", sep="_")
- sink_term <- paste("-", k_compound_sink, "*", new_boxes[[1]])
- # Only add sink term to first (or only) box for SFO and SFORB models
+ # Add first-order sink term to first (or only) box for SFO and SFORB models
if(spec[[varname]]$type %in% c("SFO", "SFORB")) {
+ k_compound_sink <- paste("k", new_boxes[[1]], "sink", sep="_")
+ sink_term <- paste("-", k_compound_sink, "*", new_boxes[[1]])
new_diffs[[1]] <- paste(new_diffs[[1]], sink_term)
+ new_parms <- k_compound_sink
+ }
+ if(spec[[varname]]$type == "FOMC") {
+ if(match(varname, obs_vars) != 1) {
+ stop("Type FOMC is only allowed for the first compartment, which is assumed to be the source compartment")
+ }
+ # From p. 53 of the FOCUS kinetics report
+ fomc_term <- paste("(alpha/beta) * ((time/beta) + 1)^-1 *", new_boxes[[1]])
+ new_diffs[[1]] <- paste(new_diffs[[1]], "-", fomc_term)
+ new_parms <- c("alpha", "beta")
}
- new_parms <- k_compound_sink
}
# Add reversible binding if appropriate
@@ -53,24 +68,41 @@ mkinmod <- function(spec = list(parent = list(type = "SFO", to = NA, sink = TRUE
}
# Transfer between compartments
- for (varname in names(spec)) {
+ for (varname in obs_vars) {
to <- spec[[varname]]$to
- if(!is.na(to)) {
+ if(!is.null(to)) {
origin_box <- switch(spec[[varname]]$type,
SFO = varname,
+ FOMC = varname,
SFORB = paste(varname, "free", sep="_"))
+ fraction_left <- NULL
for (target in to) {
target_box <- switch(spec[[target]]$type,
SFO = target,
SFORB = paste(target, "free", sep="_"))
+ if(spec[[varname]]$type %in% c("SFO", "SFORB")) {
+ k_from_to <- paste("k", origin_box, target_box, sep="_")
+ diffs[[origin_box]] <- paste(diffs[[origin_box]], "-", k_from_to, "*", origin_box)
+ diffs[[target_box]] <- paste(diffs[[target_box]], "+", k_from_to, "*", origin_box)
+ parms <- c(parms, k_from_to)
+ }
+ if(spec[[varname]]$type == "FOMC") {
+ fraction_to_target = paste("f_to", target, sep="_")
+ fraction_not_to_target = paste("(1 - ", fraction_to_target, ")", sep="")
+ if(is.null(fraction_left)) {
+ fraction_really_to_target = fraction_to_target
+ fraction_left = fraction_not_to_target
+ } else {
+ fraction_really_to_target = paste(fraction_left, " * ", fraction_to_target, sep="")
+ fraction_left = paste(fraction_left, " * ", fraction_not_to_target, sep="")
+ }
+ diffs[[target_box]] <- paste(diffs[[target_box]], "+", fraction_really_to_target, "*", fomc_term)
+ parms <- c(parms, fraction_to_target)
+ }
}
- k_from_to <- paste("k", origin_box, target_box, sep="_")
- diffs[[origin_box]] <- paste(diffs[[origin_box]], "-", k_from_to, "*", origin_box)
- diffs[[target_box]] <- paste(diffs[[target_box]], "+", k_from_to, "*", origin_box)
- parms <- c(parms, k_from_to)
}
}
model <- list(diffs = diffs, parms = parms, map = map)
class(model) <- "mkinmod"
- return(model)
+ invisible(model)
}

Contact - Imprint