diff options
Diffstat (limited to 'CakeHelpers.R')
| -rwxr-xr-x | CakeHelpers.R | 118 |
1 files changed, 118 insertions, 0 deletions
diff --git a/CakeHelpers.R b/CakeHelpers.R new file mode 100755 index 0000000..5af51cc --- /dev/null +++ b/CakeHelpers.R @@ -0,0 +1,118 @@ +# $Id$ +# Some of the CAKE R modules are based on mkin, +# Developed by Tessella Ltd for Syngenta: Copyright (C) 2011-2016 Syngenta +# Tessella Project Reference: 6245, 7247, 8361, 7414 + +# The CAKE R modules are free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more 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/>. + +# Shifts parameters slightly away from boundaries specified in "lower" and +# "upper" (to avoid computational issues after parameter transforms in modFit). +ShiftAwayFromBoundaries <- function(parameters, lower, upper) { + parametersOnLowerBound = which(parameters == lower) + parameters[parametersOnLowerBound] <- parameters[parametersOnLowerBound] * (1 + .Machine$double.eps) + .Machine$double.xmin + + parametersOnUpperBound = which(parameters == upper) + parameters[parametersOnUpperBound] <- parameters[parametersOnUpperBound] * (1 - .Machine$double.neg.eps) - .Machine$double.xmin + + return(parameters) +} + +# Adjusts stated initial values to put into the ODE solver. +# +# odeini: The initial values to adjust (in the form that would be fed into the ode function). +# cake.model: The expression of the model that we are solving. +# odeparms: The parameters for the ODE (in the form that would be fed into the ode function). +# +# Returns: Adjusted initial values. +AdjustOdeInitialValues <- function(odeini, cake.model, odeparms) { + odeini.names <- names(odeini) + + for (ini.name in odeini.names) { + # For DFOP metabolites in two compartments, need to calculate some initial conditions for the ODEs. + if (!(ini.name %in% names(cake.model$diffs))){ + subcompartment1.name <- paste(ini.name, "1", sep="_") + subcompartment2.name <- paste(ini.name, "2", sep="_") + + if (subcompartment1.name %in% names(cake.model$diffs) && subcompartment2.name %in% names(cake.model$diffs)){ + g.parameter.name = paste("g", ini.name, sep="_") + + odeini[[subcompartment1.name]] <- odeini[[ini.name]] * odeparms[[g.parameter.name]] + odeini[[subcompartment2.name]] <- odeini[[ini.name]] * (1 - odeparms[[g.parameter.name]]) + } + } + } + + # It is important that these parameters are stated in the same order as the differential equations. + return(odeini[names(cake.model$diffs)]) +} + +# Post-processes the output from the ODE solver (or analytical process), including recombination of sub-compartments. +# +# odeoutput: The output of the ODE solver. +# cake.model: The expression of the model that we are solving. +# atol: The tolerance to which the solution has been calculated. +# +# Returns: Post-processed/transformed ODE output. +PostProcessOdeOutput <- function(odeoutput, cake.model, atol) { + out_transformed <- data.frame(time = odeoutput[, "time"]) + + # Replace values that are incalculably small with 0. + for (col.name in colnames(odeoutput)) { + if (col.name == "time") { + next + } + + # If we have non-NaN, positive outputs... + if (length(odeoutput[, col.name][!is.nan(odeoutput[, col.name]) && odeoutput[, col.name] > 0]) > 0) { + # ...then replace the NaN outputs. + odeoutput[, col.name][is.nan(odeoutput[, col.name])] <- 0 + } + + # Round outputs smaller than the used tolerance down to 0. + odeoutput[, col.name][odeoutput[, col.name] < atol] <- 0 + } + + # Re-combine sub-compartments (if required) + for (compartment.name in names(cake.model$map)) { + if (length(cake.model$map[[compartment.name]]) == 1) { + out_transformed[compartment.name] <- odeoutput[, compartment.name] + } else { + out_transformed[compartment.name] <- rowSums(odeoutput[, cake.model$map[[compartment.name]]]) + } + } + + return(out_transformed) +} + +# Reorganises data in a wide format to a long format. +# +# wide_data: The data in wide format. +# time: The name of the time variable in wide_data (default "t"). +# +# Returns: Reorganised data. +wide_to_long <- function(wide_data, time = "t") { + colnames <- names(wide_data) + + if (!(time %in% colnames)) { + stop("The data in wide format have to contain a variable named ", time, ".") + } + + vars <- subset(colnames, colnames != time) + n <- length(colnames) - 1 + long_data <- data.frame(name = rep(vars, each = length(wide_data[[time]])), + time = as.numeric(rep(wide_data[[time]], n)), value = as.numeric(unlist(wide_data[vars])), + row.names = NULL) + + return(long_data) +}
\ No newline at end of file |