diff options
| -rw-r--r-- | CakeInit.R | 2 | ||||
| -rw-r--r-- | CakeIrlsPlot.R | 25 | ||||
| -rw-r--r-- | CakeOlsPlot.R | 175 |
3 files changed, 1 insertions, 201 deletions
@@ -39,7 +39,7 @@ options(error=recover) options(warn=-1) # When running from C#, this must match the C# CAKE version -CAKE.version<-"3.2" +CAKE.version<-"3.3" # The number of data points to use to draw lines on CAKE plots. CAKE.plots.resolution<-401 diff --git a/CakeIrlsPlot.R b/CakeIrlsPlot.R deleted file mode 100644 index 5a7992f..0000000 --- a/CakeIrlsPlot.R +++ /dev/null @@ -1,25 +0,0 @@ -#$Id$ -# Generates fitted curves so the GUI can plot them -# Author: Rob Nelson (Tessella) -# Tessella Project Reference: 6245 - -# Copyright (C) 2011 Syngenta -# -# This program is 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/>.” - -CakeIrlsPlot <- function(fit, xlim = range(fit$data$time), ...) -{ - CakePlotInit(fit, xlim, ...) -} - diff --git a/CakeOlsPlot.R b/CakeOlsPlot.R deleted file mode 100644 index e92742c..0000000 --- a/CakeOlsPlot.R +++ /dev/null @@ -1,175 +0,0 @@ -#$Id$ -# Generates fitted curves so the GUI can plot them -# Based on code in IRLSkinfit -# Author: Rob Nelson (Tessella) -# Modifications developed by Tessella Plc for Syngenta: Copyright (C) 2011 Syngenta -# Tessella Project Reference: 6245 -# -# This program is 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/>.” - -CakePlotInit <- function(fit, xlim = range(fit$data$time), ...) -{ - t.map.names <- names(fit$map) - metabolites <- grep("[A-Z]\\d", t.map.names, value = TRUE) - t.map.rest <- setdiff(t.map.names, metabolites) - - # Generate the normal graphs. - final <- CakeOlsPlot(fit, xlim) - final_scaled <- final - - if(length(metabolites) > 0){ - for(i in 1:length(metabolites)) - { - metabolite <- metabolites[i] - decay_var <- paste("k", metabolite, sep="_") - - # calculate the new ffm and generate the two ffm scale charts - regex <- paste("f_(.+)_to", metabolite, sep="_") - decays = grep(regex, names(fit$par), value = TRUE) - ffm_fitted <- sum(fit$par[decays]) - normal <- final - ffm_scale <- NULL - - # Generate the DT50=1000d and ffm as fitted. - k_new <- fit$par[decay_var]*fit$distimes[metabolite,]["DT50"]/1000; - fit$par[decay_var]<- k_new[metabolite,] - dt50_1000_ffm_fitted <- CakeOlsPlot(fit, xlim)[metabolite] - - naming <- c(names(final), paste(metabolite, "DT50_1000_FFM_FITTED", sep="_")) - normal <- c(final, dt50_1000_ffm_fitted) - names(normal) <- naming - final <- normal - - # Generate the scaled FFM - if(ffm_fitted != 0) - { - ffm_scale <- 1 / ffm_fitted - final_scaled <- final[c("time", metabolite, paste(metabolite, "DT50_1000_FFM_FITTED", sep="_"))] - final_scaled[t.map.rest] <- NULL; - final_frame <- as.data.frame(final_scaled) - new_names <- c(paste(metabolite, "DT50_FITTED_FFM_1", sep="_"), paste(metabolite, "DT50_1000_FFM_1", sep="_")) - names(final_frame) <- c("time", new_names) - final_frame[new_names]<-final_frame[new_names]*ffm_scale; - - cat("<PLOT MODEL START>\n") - - write.table(final_frame, quote=FALSE) - - cat("<PLOT MODEL END>\n") - } - } - } - - cat("<PLOT MODEL START>\n") - - write.table(final, quote=FALSE) - - cat("<PLOT MODEL END>\n") - - # View(final) -} - -CakeOlsPlot <- function(fit, xlim = range(fit$data$time), scale_x = 1.0, ...) -{ - solution = fit$solution - if ( is.null(solution) ) { - solution <- "deSolve" - } - atol = fit$atol - if ( is.null(atol) ) { - atol = 1.0e-6 - } - - fixed <- fit$fixed$value - names(fixed) <- rownames(fit$fixed) - parms.all <- c(fit$par, fixed) - ininames <- c( - rownames(subset(fit$start, type == "state")), - rownames(subset(fit$fixed, type == "state"))) - odeini <- parms.all[ininames] - names(odeini) <- names(fit$diffs) - - outtimes <- seq(0, xlim[2], length.out=101) * scale_x - - odenames <- c( - rownames(subset(fit$start, type == "deparm")), - rownames(subset(fit$fixed, type == "deparm"))) - odeparms <- parms.all[odenames] - - # Solve the system - evalparse <- function(string) - { - eval(parse(text=string), as.list(c(odeparms, odeini))) - } - if (solution == "analytical") { - parent.type = names(fit$map[[1]])[1] - parent.name = names(fit$diffs)[[1]] - o <- switch(parent.type, - SFO = SFO.solution(outtimes, - evalparse(parent.name), - evalparse(paste("k", parent.name, sep="_"))), - FOMC = FOMC.solution(outtimes, - evalparse(parent.name), - evalparse("alpha"), evalparse("beta")), - DFOP = DFOP.solution(outtimes, - evalparse(parent.name), - evalparse("k1"), evalparse("k2"), - evalparse("g")), - HS = HS.solution(outtimes, - evalparse(parent.name), - evalparse("k1"), evalparse("k2"), - evalparse("tb")), - SFORB = SFORB.solution(outtimes, - evalparse(parent.name), - evalparse(paste("k", parent.name, "free_bound", sep="_")), - evalparse(paste("k", parent.name, "bound_free", sep="_")), - evalparse(paste("k", parent.name, "free_sink", sep="_"))) - ) - out <- cbind(outtimes, o) - dimnames(out) <- list(outtimes, c("time", parent.name)) - } - if (solution == "eigen") { - coefmat.num <- matrix(sapply(as.vector(fit$coefmat), evalparse), - nrow = length(odeini)) - e <- eigen(coefmat.num) - c <- solve(e$vectors, odeini) - f.out <- function(t) { - e$vectors %*% diag(exp(e$values * t), nrow=length(odeini)) %*% c - } - o <- matrix(mapply(f.out, outtimes), - nrow = length(odeini), ncol = length(outtimes)) - dimnames(o) <- list(names(odeini), NULL) - out <- cbind(time = outtimes, t(o)) - } - if (solution == "deSolve") { - out <- ode( - y = odeini, - times = outtimes, - func = fit$mkindiff, - parms = odeparms, - atol = atol - ) - } - - # Output transformation for models with unobserved compartments like SFORB - out_transformed <- data.frame(time = out[,"time"]) - for (var in names(fit$map)) { - if(length(fit$map[[var]]) == 1) { - out_transformed[var] <- out[, var] - } else { - out_transformed[var] <- rowSums(out[, fit$map[[var]]]) - } - } - return(out_transformed) -} |