diff options
| -rw-r--r-- | R/IORE.solution.R | 4 | ||||
| -rw-r--r-- | R/mkinfit.R | 4 | ||||
| -rw-r--r-- | R/mkinmod.R | 28 | ||||
| -rw-r--r-- | R/mkinpredict.R | 6 | ||||
| -rw-r--r-- | R/transform_odeparms.R | 14 | ||||
| -rw-r--r-- | man/FOMC.solution.Rd | 2 | ||||
| -rw-r--r-- | man/IORE.solution.Rd | 34 | ||||
| -rw-r--r-- | man/mkinmod.Rd | 5 |
8 files changed, 86 insertions, 11 deletions
diff --git a/R/IORE.solution.R b/R/IORE.solution.R new file mode 100644 index 00000000..9546ce56 --- /dev/null +++ b/R/IORE.solution.R @@ -0,0 +1,4 @@ +IORE.solution <- function(t, parent.0, k.iore, N)
+{
+ parent = (parent.0^(1 - N) - (1 - N) * k.iore * t)^(1/(1 - N))
+}
diff --git a/R/mkinfit.R b/R/mkinfit.R index b7ca1d74..26fc0e6b 100644 --- a/R/mkinfit.R +++ b/R/mkinfit.R @@ -90,13 +90,15 @@ mkinfit <- function(mkinmod, observed, defaultpar.names <- setdiff(mkinmod$parms, names(parms.ini))
for (parmname in defaultpar.names) {
# Default values for rate constants, depending on the parameterisation
- if (substr(parmname, 1, 2) == "k_") {
+ if (grepl("^k", parmname)) {
parms.ini[parmname] = 0.1 + k_salt
k_salt = k_salt + 1e-4
}
# Default values for rate constants for reversible binding
if (grepl("free_bound$", parmname)) parms.ini[parmname] = 0.1
if (grepl("bound_free$", parmname)) parms.ini[parmname] = 0.02
+ # Default values for IORE exponents
+ if (grepl("^N", parmname)) parms.ini[parmname] = 1
# Default values for the FOMC, DFOP and HS models
if (parmname == "alpha") parms.ini[parmname] = 1
if (parmname == "beta") parms.ini[parmname] = 10
diff --git a/R/mkinmod.R b/R/mkinmod.R index eb290719..321887fc 100644 --- a/R/mkinmod.R +++ b/R/mkinmod.R @@ -41,9 +41,8 @@ mkinmod <- function(..., use_of_ff = "min", speclist = NULL) map <- list()
# }}}
- # Do not return a coefficient matrix mat when FOMC, DFOP or HS are used for
- # the parent compound {{{
- if(spec[[1]]$type %in% c("FOMC", "DFOP", "HS")) {
+ # Do not return a coefficient matrix mat when FOMC, IORE, DFOP or HS is used for the parent {{{
+ if(spec[[1]]$type %in% c("FOMC", "IORE", "DFOP", "HS")) {
mat = FALSE
} else mat = TRUE
#}}}
@@ -55,8 +54,8 @@ mkinmod <- function(..., use_of_ff = "min", speclist = NULL) # Check the type component of the compartment specification {{{
if(is.null(spec[[varname]]$type)) stop(
"Every part of the model specification must be a list containing a type component")
- if(!spec[[varname]]$type %in% c("SFO", "FOMC", "DFOP", "HS", "SFORB")) stop(
- "Available types are SFO, FOMC, DFOP, HS and SFORB only")
+ if(!spec[[varname]]$type %in% c("SFO", "FOMC", "IORE", "DFOP", "HS", "SFORB")) stop(
+ "Available types are SFO, FOMC, IORE, DFOP, HS and SFORB only")
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"))
@@ -65,6 +64,7 @@ mkinmod <- function(..., use_of_ff = "min", speclist = NULL) new_boxes <- switch(spec[[varname]]$type,
SFO = varname,
FOMC = varname,
+ IORE = varname,
DFOP = varname,
HS = varname,
SFORB = paste(varname, c("free", "bound"), sep = "_")
@@ -85,20 +85,36 @@ mkinmod <- function(..., use_of_ff = "min", speclist = NULL) # Turn on sink if this is not explicitly excluded by the user by
# specifying sink=FALSE
if(is.null(spec[[varname]]$sink)) spec[[varname]]$sink <- TRUE
- if(spec[[varname]]$type %in% c("SFO", "SFORB")) { # {{{ Add SFO or SFORB decline term
+ if(spec[[varname]]$type %in% c("SFO", "IORE", "SFORB")) { # {{{ Add decline term
if (use_of_ff == "min") { # Minimum use of formation fractions
if(spec[[varname]]$sink) {
# If sink is required, add first-order sink term
k_compound_sink <- paste("k", box_1, "sink", sep = "_")
+ if(spec[[varname]]$type == "IORE") {
+ k_compound_sink <- paste("k.iore", box_1, "sink", sep = "_")
+ }
parms <- c(parms, k_compound_sink)
decline_term <- paste(k_compound_sink, "*", box_1)
+ if(spec[[varname]]$type == "IORE") {
+ N <- paste("N", box_1, sep = "_")
+ parms <- c(parms, N)
+ decline_term <- paste0(decline_term, "^", N)
+ }
} else { # otherwise no decline term needed here
decline_term = "0"
}
} else {
k_compound <- paste("k", box_1, sep = "_")
+ if(spec[[varname]]$type == "IORE") {
+ k_compound <- paste("k.iore", box_1, sep = "_")
+ }
parms <- c(parms, k_compound)
decline_term <- paste(k_compound, "*", box_1)
+ if(spec[[varname]]$type == "IORE") {
+ N <- paste("N", box_1, sep = "_")
+ parms <- c(parms, N)
+ decline_term <- paste0(decline_term, "^", N)
+ }
}
} #}}}
if(spec[[varname]]$type == "FOMC") { # {{{ Add FOMC decline term
diff --git a/R/mkinpredict.R b/R/mkinpredict.R index 14efc568..43b5d2f3 100644 --- a/R/mkinpredict.R +++ b/R/mkinpredict.R @@ -50,6 +50,12 @@ mkinpredict <- function(mkinmod, odeparms, odeini, FOMC = FOMC.solution(outtimes, evalparse(parent.name), evalparse("alpha"), evalparse("beta")), + IORE = IORE.solution(outtimes, + evalparse(parent.name), + ifelse(mkinmod$use_of_ff == "min", + evalparse(paste("k.iore", parent.name, "sink", sep="_")), + evalparse(paste("k.iore", parent.name, sep="_"))), + evalparse("N_parent")), DFOP = DFOP.solution(outtimes, evalparse(parent.name), evalparse("k1"), evalparse("k2"), diff --git a/R/transform_odeparms.R b/R/transform_odeparms.R index 912a5c0a..a36b7eae 100644 --- a/R/transform_odeparms.R +++ b/R/transform_odeparms.R @@ -33,6 +33,8 @@ transform_odeparms <- function(parms, mkinmod, # Log transformation for rate constants if requested
k <- parms[grep("^k_", names(parms))]
+ k.iore <- parms[grep("^k.iore_", names(parms))]
+ k <- c(k, k.iore)
if (length(k) > 0) {
if(transform_rates) {
transparms[paste0("log_", names(k))] <- log(k)
@@ -40,6 +42,10 @@ transform_odeparms <- function(parms, mkinmod, else transparms[names(k)] <- k
}
+ # Do not transform exponents in IORE models
+ N <- parms[grep("^N", names(parms))]
+ transparms[names(N)] <- N
+
# Go through state variables and apply isometric logratio transformation to
# formation fractions if requested
mod_vars = names(spec)
@@ -98,8 +104,10 @@ backtransform_odeparms <- function(transparms, mkinmod, # Exponential transformation for rate constants
if(transform_rates) {
trans_k <- transparms[grep("^log_k_", names(transparms))]
+ trans_k.iore <- transparms[grep("^log_k.iore_", names(transparms))]
+ trans_k = c(trans_k, trans_k.iore)
if (length(trans_k) > 0) {
- k_names <- gsub("^log_k_", "k_", names(trans_k))
+ k_names <- gsub("^log_k", "k", names(trans_k))
parms[k_names] <- exp(trans_k)
}
} else {
@@ -107,6 +115,10 @@ backtransform_odeparms <- function(transparms, mkinmod, parms[names(trans_k)] <- trans_k
}
+ # Do not transform exponents in IORE models
+ N <- transparms[grep("^N", names(transparms))]
+ parms[names(N)] <- N
+
# Go through state variables and apply inverse isometric logratio transformation
mod_vars = names(spec)
for (box in mod_vars) {
diff --git a/man/FOMC.solution.Rd b/man/FOMC.solution.Rd index d04d34e1..f4a26e41 100644 --- a/man/FOMC.solution.Rd +++ b/man/FOMC.solution.Rd @@ -44,6 +44,6 @@ FOMC.solution(t, parent.0, alpha, beta) Technology} \bold{24}, 1032-1038
}
\examples{
- \dontrun{plot(function(x) FOMC.solution(x, 100, 10, 2), 0, 2)}
+ \dontrun{plot(function(x) FOMC.solution(x, 100, 10, 2), 0, 2, ylim = c(0, 100))}
}
\keyword{ manip }
diff --git a/man/IORE.solution.Rd b/man/IORE.solution.Rd new file mode 100644 index 00000000..65dac9a7 --- /dev/null +++ b/man/IORE.solution.Rd @@ -0,0 +1,34 @@ +\name{IORE.solution}
+\Rdversion{1.1}
+\alias{IORE.solution}
+\title{ Indeterminate order rate equation kinetics }
+\description{
+ Function describing exponential decline from a defined starting value, with
+ a concentration dependent rate constant.
+}
+\usage{
+ IORE.solution(t, parent.0, k.iore, N)
+}
+\arguments{
+ \item{t}{ Time. }
+ \item{parent.0}{ Starting value for the response variable at time zero. }
+ \item{k.iore}{ Rate constant. Note that this depends on the concentration units used. }
+ \item{N}{ Exponent describing the nonlinearity of the rate equation }
+}
+\note{
+ The solution of the IORE kinetic model reduces to the
+ \code{\link{SFO.solution}} if N = 1.
+}
+\value{
+ The value of the response variable at time \code{t}.
+}
+\references{
+ NAFTA Technical Working Group on Pesticides (not dated) Guidance for
+ Evaluating and Calculating Degradation Kinetics in Environmental
+ Media
+}
+\examples{
+ \dontrun{plot(function(x) IORE.solution(x, 100, 0.2, 1.3), 0, 2,
+ ylim = c(0, 100))}
+}
+\keyword{ manip }
diff --git a/man/mkinmod.Rd b/man/mkinmod.Rd index 76127c58..63087d49 100644 --- a/man/mkinmod.Rd +++ b/man/mkinmod.Rd @@ -18,7 +18,7 @@ mkinmod(..., use_of_ff = "min", speclist = NULL) for the variable. Currently, single first order kinetics "SFO" or
single first order with reversible binding "SFORB" are implemented for all
variables, while
- "FOMC", "DFOP" and "HS" can additionally be chosen for the first
+ "FOMC", "IORE", "DFOP" and "HS" can additionally be chosen for the first
variable which is assumed to be the source compartment.
Additionally, each component of the list can include a character vector \code{to},
specifying names of variables to which a transfer is to be assumed in the
@@ -46,7 +46,8 @@ mkinmod(..., use_of_ff = "min", speclist = NULL) \item{map}{ A list containing named character vectors for each observed variable, specifying
the modelling variables by which it is represented. }
\item{use_of_ff}{ The content of \code{use_of_ff} is passed on in this list component. }
- \item{coefmat}{ The coefficient matrix, if the system of differential equations can be represented by one. }
+ \item{coefmat}{ The coefficient matrix, if the system of differential equations can be
+ represented by one. }
}
\author{
Johannes Ranke
|