From 3a6b9f52c74d6ef88a8d32c50e42864b3f251719 Mon Sep 17 00:00:00 2001 From: jranke Date: Thu, 15 Mar 2012 15:54:14 +0000 Subject: Update kinfit and mkin to the latest version published on BerliOS. git-svn-id: svn+ssh://svn.r-forge.r-project.org/svnroot/kinfit/pkg/kinfit@17 edb9625f-4e0d-4859-8d74-9fd3b1da38cb --- R/kinerrmin.R | 20 ++++++ R/kinfit.R | 115 +++++++++++++++++++++++++++------- R/kinobject.R | 20 ++++++ R/kinobjects.R | 20 ++++++ R/kinplot.R | 128 ++++++++++++++++++++++++++------------ R/kinreport.R | 97 +++++++++++++++++++++-------- R/kinresplot.R | 23 ++++++- R/kinresults.R | 175 ++++++++++++++++++++++++++++++++++++---------------- R/kinwrite.KinGUI.R | 20 ++++++ 9 files changed, 478 insertions(+), 140 deletions(-) (limited to 'R') diff --git a/R/kinerrmin.R b/R/kinerrmin.R index 8ce30a3..dbd64b3 100644 --- a/R/kinerrmin.R +++ b/R/kinerrmin.R @@ -1,3 +1,23 @@ +# $Id: kinerrmin.R 59 2010-07-28 12:29:15Z jranke $ + +# Copyright (C) 2008-2010 Johannes Ranke +# Contact: mkin-devel@lists.berlios.de + +# This file is part of the R package kinfit + +# kinfit 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 + kinerrmin <- function(kinfits, kinmodel = "SFO", alpha = 0.05) { m = kinfits[[kinmodel]] diff --git a/R/kinfit.R b/R/kinfit.R index 42a6520..d55180e 100644 --- a/R/kinfit.R +++ b/R/kinfit.R @@ -1,43 +1,76 @@ +# $Id: kinfit.R 116 2011-06-14 08:46:47Z kati $ + +# Copyright (C) 2008-2010 Johannes Ranke +# Contact: mkin-devel@lists.berlios.de + +# This file is part of the R package kinfit + +# kinfit 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 + kinfit <- function(kindata, kinmodels = c("SFO"), parent.0.user = NA, + parent.0.fixed = FALSE, start.SFO = list(parent.0 = NA, k = NA), start.FOMC = list(parent.0 = NA, alpha = NA, beta = NA), start.DFOP = list(parent.0 = NA, k1 = NA, k2 = NA, g = NA), start.HS = list(parent.0 = NA, k1 = NA, k2 = NA, tb = NA), - algorithm = "port") + algorithm = "default") { + kindata <- subset(kindata, !is.na(kindata$parent)) kinfits <- list() if (!is.na(parent.0.user)) { start.SFO$parent.0 = parent.0.user start.FOMC$parent.0 = parent.0.user + start.DFOP$parent.0 = parent.0.user + start.HS$parent.0 = parent.0.user } lmlogged = lm(log(parent) ~ t, data = kindata) + k.est = -coef(lmlogged)[["t"]] for (kinmodel in kinmodels) { if (kinmodel == "SFO") { if (is.na(start.SFO$parent.0)) { - start.SFO$parent.0 = max(kindata$parent) + start.SFO$parent.0 = max(kindata$parent) } if (is.na(start.SFO$k)) { start.SFO$k = - coef(lmlogged)[["t"]] } - kinfits[[kinmodel]] = try( - nls(parent ~ SFO(t, parent.0, k), - data = kindata, model = TRUE, - start = start.SFO, - algorithm = algorithm), silent=TRUE) + if (parent.0.fixed) + { + start.SFO = start.SFO[-1] + kinfits[[kinmodel]] = try( + nls(parent ~ SFO(t, parent.0.user, k), + data = kindata, model = TRUE, + start = start.SFO, + algorithm = algorithm), silent=TRUE) + } else { + kinfits[[kinmodel]] = try( + nls(parent ~ SFO(t, parent.0, k), + data = kindata, model = TRUE, + start = start.SFO, + algorithm = algorithm), silent=TRUE) + } + k.est = coef(kinfits$SFO)[["k"]] } - k.est = ifelse(is.na(coef(kinfits$SFO)[["k"]]), - -coef(lmlogged)[["t"]], - coef(kinfits$SFO)[["k"]]) if (kinmodel == "FOMC") { if (is.na(start.FOMC$parent.0)) { - start.FOMC$parent.0 = max(kindata$parent) + start.FOMC$parent.0 = max(kindata$parent) } if (is.na(start.FOMC$alpha)) { start.FOMC$alpha = 1 @@ -45,15 +78,29 @@ kinfit <- function(kindata, kinmodels = c("SFO"), if (is.na(start.FOMC$beta)) { start.FOMC$beta = start.FOMC$alpha / k.est } - kinfits[[kinmodel]] = try( - nls(parent ~ FOMC(t, parent.0, alpha, beta), - data = kindata, model = TRUE, - start = start.FOMC, - algorithm = algorithm), silent=TRUE) + + if (parent.0.fixed) + { + start.FOMC = list(alpha = start.FOMC$alpha, beta = start.FOMC$beta) + + kinfits[[kinmodel]] = try( + nls(parent ~ FOMC(t, parent.0.user, alpha, beta), + data = kindata, model = TRUE, + start = start.FOMC, + algorithm = algorithm), silent=TRUE) + + } else { + kinfits[[kinmodel]] = try( + nls(parent ~ FOMC(t, parent.0, alpha, beta), + data = kindata, model = TRUE, + start = start.FOMC, + algorithm = algorithm), silent=TRUE) + + } } if (kinmodel == "DFOP") { if (is.na(start.DFOP$parent.0)) { - start.DFOP$parent.0 = max(kindata$parent) + start.DFOP$parent.0 = max(kindata$parent) } if (is.na(start.DFOP$k1)) { start.DFOP$k1 = k.est * 2 @@ -64,15 +111,26 @@ kinfit <- function(kindata, kinmodels = c("SFO"), if (is.na(start.DFOP$g)) { start.DFOP$g = 0.5 } - kinfits[[kinmodel]] = try( + if (parent.0.fixed) + { + start.DFOP = list(k1 = start.DFOP$k1, k2 = start.DFOP$k2, g = start.DFOP$g) + + kinfits[[kinmodel]] = try( + nls(parent ~ DFOP(t, parent.0.user, k1, k2, g), + data = kindata, model = TRUE, + start = start.DFOP, + algorithm = algorithm), silent=TRUE) + }else{ + kinfits[[kinmodel]] = try( nls(parent ~ DFOP(t, parent.0, k1, k2, g), data = kindata, model = TRUE, start = start.DFOP, - algorithm = algorithm), silent=TRUE) + algorithm = algorithm), silent=TRUE) + } } if (kinmodel == "HS") { if (is.na(start.HS$parent.0)) { - start.HS$parent.0 = max(kindata$parent) + start.HS$parent.0 = max(kindata$parent) } if (is.na(start.HS$k1)) { start.HS$k1 = k.est @@ -83,11 +141,24 @@ kinfit <- function(kindata, kinmodels = c("SFO"), if (is.na(start.HS$tb)) { start.HS$tb = 0.05 * max(kindata$t) } - kinfits[[kinmodel]] = try( + + if (parent.0.fixed) + { + + start.HS = list(k1 = start.HS$k1, k2 = start.HS$k2, tb = start.HS$tb) + + kinfits[[kinmodel]] = try( + nls(parent ~ HS(t, parent.0.user, k1, k2, tb), + data = kindata, model = TRUE, + start = start.HS, + algorithm = algorithm), silent=TRUE) + }else{ + kinfits[[kinmodel]] = try( nls(parent ~ HS(t, parent.0, k1, k2, tb), data = kindata, model = TRUE, start = start.HS, - algorithm = algorithm), silent=TRUE) + algorithm = algorithm), silent=TRUE) + } } } return(kinfits) diff --git a/R/kinobject.R b/R/kinobject.R index de6f6af..2ecc962 100644 --- a/R/kinobject.R +++ b/R/kinobject.R @@ -1,3 +1,23 @@ +# $Id: kinobject.R 59 2010-07-28 12:29:15Z jranke $ + +# Copyright (C) 2008-2010 Johannes Ranke +# Contact: mkin-devel@lists.berlios.de + +# This file is part of the R package kinfit + +# kinfit 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 + kinobject <- function(parent, type, system, layers = NA, sampling_times = NA) { diff --git a/R/kinobjects.R b/R/kinobjects.R index a767cea..33abe61 100644 --- a/R/kinobjects.R +++ b/R/kinobjects.R @@ -1,3 +1,23 @@ +# $Id: kinobjects.R 59 2010-07-28 12:29:15Z jranke $ + +# Copyright (C) 2008-2010 Johannes Ranke +# Contact: mkin-devel@lists.berlios.de + +# This file is part of the R package kinfit + +# kinfit 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 + kinobjects<- function(parent, type, systems, layers = NA, sampling_times = NA) { diff --git a/R/kinplot.R b/R/kinplot.R index ace1270..394e271 100644 --- a/R/kinplot.R +++ b/R/kinplot.R @@ -1,61 +1,113 @@ +# $Id: kinplot.R 117 2011-06-14 08:52:14Z kati $ + +# Copyright (C) 2008-2010 Johannes Ranke +# Contact: mkin-devel@lists.berlios.de + +# This file is part of the R package kinfit + +# kinfit 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 + kinplot <- function(kinobject, + main = "", xlab = "Time [days]", ylab = "Parent [% of applied radioactivity]", ylim = c("auto", "auto"), lpos = "topright") { kindata <- na.omit(kinobject$data) kinfits <- kinobject$fits - if (ylim[1] == "auto") ylim[1] <- 0 - if (ylim[2] == "auto") ylim[2] <- max(kindata$parent) - ylim <- as.numeric(ylim) + if (ylim[1] == "auto") ylim[1] <- 0 + if (ylim[2] == "auto") ylim[2] <- max(kindata$parent) + ylim <- as.numeric(ylim) plot(kindata$t, kindata$parent, + main = main, xlab = xlab, ylab = ylab, ylim = ylim - ) + ) n.m <- length(kinfits) colors <- ltys <- 1:n.m names(colors) <- names(ltys) <- names(kinfits) - ltext <- paste(kinobject$parent, "measured") + ltext <- paste(kinobject$parent, "measured") for (kinmodel in names(kinfits)) { m = kinfits[[kinmodel]] if(class(m) == "nls") { - switch(kinmodel, - SFO = curve(SFO(x, - coef(m)[["parent.0"]], - coef(m)[["k"]]), - from = min(kindata$t), to = max(kindata$t), add=TRUE, - col = colors[[kinmodel]], - lty = ltys[[kinmodel]]), - FOMC = curve(FOMC(x, - coef(m)[["parent.0"]], - coef(m)[["alpha"]], - coef(m)[["beta"]]), - from = min(kindata$t), to = max(kindata$t), add=TRUE, - col = colors[[kinmodel]], - lty = ltys[[kinmodel]]), - HS = curve(HS(x, - coef(m)[["parent.0"]], - coef(m)[["k1"]], - coef(m)[["k2"]], - coef(m)[["tb"]]), - from = min(kindata$t), to = max(kindata$t), add=TRUE, - col = colors[[kinmodel]], - lty = ltys[[kinmodel]]), - DFOP = curve(DFOP(x, - coef(m)[["parent.0"]], - coef(m)[["k1"]], - coef(m)[["k2"]], - coef(m)[["g"]]), - from = min(kindata$t), to = max(kindata$t), add=TRUE, - col = colors[[kinmodel]], - lty = ltys[[kinmodel]])) - ltext <- c(ltext, paste("Fitted", kinmodel, "model")) + if (!"parent.0" %in% names(coef(m))) { + switch(kinmodel, + SFO = lines( + t <- seq(min(kindata$t), max(kindata$t), length.out=500), + predict(m, + newdata = data.frame(t)), + col = colors[[kinmodel]], + lty = ltys[[kinmodel]]), + FOMC = lines( + t <- seq(min(kindata$t), max(kindata$t), length.out=500), + predict(m, + newdata = data.frame(t)), + col = colors[[kinmodel]], + lty = ltys[[kinmodel]]), + HS = lines( + t <- seq(min(kindata$t), max(kindata$t), length.out=500), + predict(m, + newdata = data.frame(t)), + col = colors[[kinmodel]], + lty = ltys[[kinmodel]]), + DFOP = lines( + t <- seq(min(kindata$t), max(kindata$t), length.out=500), + predict(m, + newdata = data.frame(t)), + col = colors[[kinmodel]], + lty = ltys[[kinmodel]]) + ) + ltext <- c(ltext, paste("Fitted", kinmodel, "model")) + } else { + switch(kinmodel, + SFO = curve(SFO(x, + coef(m)[["parent.0"]], + coef(m)[["k"]]), + from = min(kindata$t), to = max(kindata$t), add=TRUE, + col = colors[[kinmodel]], + lty = ltys[[kinmodel]]), + FOMC = curve(FOMC(x, + coef(m)[["parent.0"]], + coef(m)[["alpha"]], + coef(m)[["beta"]]), + from = min(kindata$t), to = max(kindata$t), add=TRUE, + col = colors[[kinmodel]], + lty = ltys[[kinmodel]]), + HS = curve(HS(x, + coef(m)[["parent.0"]], + coef(m)[["k1"]], + coef(m)[["k2"]], + coef(m)[["tb"]]), + from = min(kindata$t), to = max(kindata$t), add=TRUE, + col = colors[[kinmodel]], + lty = ltys[[kinmodel]]), + DFOP = curve(DFOP(x, + coef(m)[["parent.0"]], + coef(m)[["k1"]], + coef(m)[["k2"]], + coef(m)[["g"]]), + from = min(kindata$t), to = max(kindata$t), add=TRUE, + col = colors[[kinmodel]], + lty = ltys[[kinmodel]])) + ltext <- c(ltext, paste("Fitted", kinmodel, "model")) + } } else { - ltext <- c(ltext, paste(kinmodel, "model failed")) - ltys[[kinmodel]] <- NA + ltext <- c(ltext, paste(kinmodel, "model failed")) + ltys[[kinmodel]] <- NA } } legend(lpos, bty="n", inset = 0.05, diff --git a/R/kinreport.R b/R/kinreport.R index 4156803..d405b8e 100644 --- a/R/kinreport.R +++ b/R/kinreport.R @@ -1,43 +1,86 @@ -kinreport <- function(kinobject, file = NA, vcov = FALSE, endpoint.digits = 1) +# $Id: kinreport.R 123 2011-11-01 12:26:41Z jranke $ + +# Copyright (C) 2008-2010 Johannes Ranke +# Contact: mkin-devel@lists.berlios.de + +# This file is part of the R package kinfit + +# kinfit 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 + +kinreport <- function(kinobject, file = NA, data = TRUE, R2 = FALSE, vcov = FALSE, endpoint.digits = 1) { if (!is.na(file)) { sink(file, split=TRUE) } cat("Parent compound: ", kinobject$parent, "\n") - if (!is.null(kinobject$label)) cat("Label position:\t\t", kinobject$label, "\n") + if (!is.null(kinobject$label)) { + cat("Label position:\t\t", kinobject$label, "\n") + } cat("Study type: ", kinobject$type, "\n") cat("System: ", kinobject$system, "\n") if (!is.null(kinobject$source)) { - cat("Source: ", kinobject$source, "\n") - } + cat("Source: ", kinobject$source, "\n") + } + cat("kinfit version: ", as.character(packageVersion("kinfit")), "\n") + cat("R version: ", paste(R.version$major, R.version$minor, sep="."), "\n") + cat("Report generated:", date(), "\n") cat("\n") + if (data) { + cat("Data:\n") + print(kinobject$data) + cat("\n") + } fit.names <- names(kinobject$fits) for (kinmodel in fit.names) { - m <- kinobject$fits[[kinmodel]] - if (!(class(m) == "try-error")) { - cat("\n\n---\n") - cat("Nonlinear least squares fit of the", kinmodel, "model\n\n") - cat("Parameter estimation:\t") - s <- summary(m) - df <- s$df[2] - p <- 1 - pt(s$parameters[,3], df = df) - parms <- cbind(s$parameters[,c(1,2,3)], "Pr(>t)" = p) - cat("\n") - print(parms, digits=3) - cat("\n") - if(vcov) - { - cat("Variance-covariance matrix:\n") - print(vcov(m)) - cat("\n") - } - cat("Chi2 error estimation:\t", - round(100 * kinobject$results$stats[kinmodel, "err.min"], digits=2), - " %\n", sep="") - cat("\n") - } + m <- kinobject$fits[[kinmodel]] + if (!(class(m) == "try-error")) { + cat("\n\n---\n") + cat("Nonlinear least squares fit of the", kinmodel, "model\n") + if (!"parent.0" %in% names(coef(m))) { + cat(paste("Initial value of parent fixed to ", m$model$parent.0.user, "\n", sep="")) + } + cat("\n") + cat("Parameter estimation:\t") + s <- summary(m) + df <- s$df[2] + p <- 1 - pt(s$parameters[,3], df = df) + parms <- matrix(nrow = nrow(s$parameters), ncol=4) + dimnames(parms) = list(rownames(s$parameters), + c("Estimate", "Std. Error", "t value", "Pr(>t)")) + parms[, c(1,2,3)] <- s$parameters[,c(1,2,3)] + parms[, 4] <- p + cat("\n") + print(parms, digits=3) + cat("\n") + if(vcov) + { + cat("Variance-covariance matrix:\n") + print(vcov(m)) + cat("\n") + } + cat("Chi2 error estimation: ", + round(100 * kinobject$results$stats[kinmodel, "err.min"], digits=2), + " %\n", sep="") + cat("\n") + if(R2) + { + cat("Coefficient of determination R2: ", + round(kinobject$results$stats[kinmodel, "R2"], digits=3), "\n") + } + } } cat("\n\n---\n") cat("Endpoint estimates\n\n") diff --git a/R/kinresplot.R b/R/kinresplot.R index be0a85d..36eb55a 100644 --- a/R/kinresplot.R +++ b/R/kinresplot.R @@ -1,3 +1,23 @@ +# $Id: kinresplot.R 106 2011-05-12 10:55:37Z jranke $ + +# Copyright (C) 2008-2010 Johannes Ranke +# Contact: mkin-devel@lists.berlios.de + +# This file is part of the R package kinfit + +# kinfit 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 + kinresplot <- function(kinobject, kinmodel, xlab = "Time [days]", ylab = "Residual [% of applied radioactivity]", maxabs = "auto") @@ -5,10 +25,11 @@ kinresplot <- function(kinobject, kinmodel, m <- kinobject$fits[[kinmodel]] t <- m$model$t residuals <- residuals(m) - if (maxabs == "auto") maxabs = max(abs(residuals)) + if (maxabs == "auto") maxabs = max(abs(residuals), na.rm = TRUE) plot(t, residuals, xlab = xlab, ylab = ylab, ylim = c( -1.2 * maxabs, 1.2 * maxabs)) + abline(h=0, lty=2) title(paste("Residuals of", kinmodel, "fit"), font.main = 1) } diff --git a/R/kinresults.R b/R/kinresults.R index 6bbff28..27bdafb 100644 --- a/R/kinresults.R +++ b/R/kinresults.R @@ -1,74 +1,145 @@ -kinresults <- function(kinfits, alpha = 0.05, SFORB=TRUE) +# $Id: kinresults.R 124 2011-11-09 10:58:27Z jranke $ + +# Copyright (C) 2008-2011 Johannes Ranke +# Contact: mkin-devel@lists.berlios.de + +# This file is part of the R package kinfit + +# kinfit 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 + +kinresults <- function(kinfits, alpha = 0.05, DTmax = 1000, SFORB=TRUE) { - kindata <- data.frame(t = kinfits[[1]]$model$t, parent = kinfits[[1]]$model$parent) - kindata.means <- aggregate(kindata, list(kindata$t), mean) + kindata <- data.frame( + t = kinfits[[1]]$model$t, + parent = kinfits[[1]]$model$parent) + kindata.means <- aggregate(kindata, list(kindata$t), mean) kindata.means.mean <- mean(kindata.means$parent, na.rm=TRUE) n.times <- length(kindata.means$parent) parms <- list() - df <- err.min <- RSS <- vector() + df <- err.min <- R2 <- RSS <- TSS <- RSS.means <- TSS.means <- vector() DT50 <- DT90 <- vector() f <- list() + for (kinmodel in names(kinfits)) { m = kinfits[[kinmodel]] if(class(m) == "nls") { kindata.means$est <- predict(m, kindata.means) - parms[[kinmodel]] <- switch(kinmodel, - SFO = list(parent.0 = coef(m)[["parent.0"]], - k = coef(m)[["k"]]), - FOMC = list(parent.0 = coef(m)[["parent.0"]], - alpha = coef(m)[["alpha"]], - beta = coef(m)[["beta"]]), - HS = list(parent.0 = coef(m)[["parent.0"]], - k1 = coef(m)[["k1"]], - k2 = coef(m)[["k2"]], - tb = coef(m)[["tb"]]), - DFOP = list(parent.0 = coef(m)[["parent.0"]], - k1 = coef(m)[["k1"]], - k2 = coef(m)[["k2"]], - g = coef(m)[["g"]])) - if(kinmodel == "DFOP" & SFORB) { - k1 = coef(m)[["k1"]] - k2 = coef(m)[["k2"]] - g = coef(m)[["g"]] - parms[["SFORB"]] = - list(parent.0 = coef(m)[["parent.0"]], - k1out = g * k1 + (1 - g) * k2, - k21 = k1 * k2 / (g * k1 + (1 - g) * k2), - k12 = (g * (1 - g) * (k1 - k2)^2) / (g * k1 + (1 - g) * k2)) + + if (!"parent.0" %in% names(coef(m))) { + + parms[[kinmodel]] <- switch(kinmodel, + SFO = list(k = coef(m)[["k"]]), + FOMC = list(alpha = coef(m)[["alpha"]], + beta = coef(m)[["beta"]]), + HS = list(k1 = coef(m)[["k1"]], + k2 = coef(m)[["k2"]], + tb = coef(m)[["tb"]]), + DFOP= list(k1 = coef(m)[["k1"]], + k2 = coef(m)[["k2"]], + g = coef(m)[["g"]])) + } else { + parms[[kinmodel]] <- switch(kinmodel, + SFO = list(parent.0 = coef(m)[["parent.0"]], + k = coef(m)[["k"]]), + FOMC = list(parent.0 = coef(m)[["parent.0"]], + alpha = coef(m)[["alpha"]], + beta = coef(m)[["beta"]]), + HS = list(parent.0 = coef(m)[["parent.0"]], + k1 = coef(m)[["k1"]], + k2 = coef(m)[["k2"]], + tb = coef(m)[["tb"]]), + DFOP = list(parent.0 = coef(m)[["parent.0"]], + k1 = coef(m)[["k1"]], + k2 = coef(m)[["k2"]], + g = coef(m)[["g"]])) + if(kinmodel == "DFOP" & SFORB) { + k1 = coef(m)[["k1"]] + k2 = coef(m)[["k2"]] + g = coef(m)[["g"]] + parms[["SFORB"]] = + list(parent.0 = coef(m)[["parent.0"]], + k1out = g * k1 + (1 - g) * k2, + k21 = k1 * k2 / (g * k1 + (1 - g) * k2), + k12 = (g * (1 - g) * (k1 - k2)^2) / (g * k1 + (1 - g) * k2)) + } + } + + n.parms = length(coef(m)) + if (!"parent.0" %in% names(coef(m))) { + + f[[kinmodel]] = switch(kinmodel, + HS = function(t, x) { + (HS(t, kinfits[[kinmodel]]$model[["parent.0.user"]], + coef(m)[["k1"]], coef(m)[["k2"]], coef(m)[["tb"]]) - + (1 - x/100) * kinfits[[kinmodel]]$model[["parent.0.user"]])^2 + }, + DFOP = function(t, x) { + (DFOP(t, kinfits[[kinmodel]]$model[["parent.0.user"]], + coef(m)[["k1"]], coef(m)[["k2"]], coef(m)[["g"]]) - + (1 - x/100) * kinfits[[kinmodel]]$model[["parent.0.user"]])^2 + } + ) + }else{ + f[[kinmodel]] = switch(kinmodel, + HS = function(t, x) { + (HS(t, coef(m)[["parent.0"]], + coef(m)[["k1"]], coef(m)[["k2"]], coef(m)[["tb"]]) - + (1 - x/100) * coef(m)[["parent.0"]])^2 + }, + DFOP = function(t, x) { + (DFOP(t, coef(m)[["parent.0"]], + coef(m)[["k1"]], coef(m)[["k2"]], coef(m)[["g"]]) - + (1 - x/100) * coef(m)[["parent.0"]])^2 + } + ) } - n.parms = length(coef(m)) - f[[kinmodel]] = switch(kinmodel, - HS = function(t, x) { - (HS(t, coef(m)[["parent.0"]], - coef(m)[["k1"]], coef(m)[["k2"]], coef(m)[["tb"]]) - - (1 - x/100) * coef(m)[["parent.0"]])^2 - }, - DFOP = function(t, x) { - (DFOP(t, coef(m)[["parent.0"]], - coef(m)[["k1"]], coef(m)[["k2"]], coef(m)[["g"]]) - - (1 - x/100) * coef(m)[["parent.0"]])^2 - } - ) + coef(m) df[[kinmodel]] = n.times - n.parms RSS[[kinmodel]] = sum(summary(m)$residuals^2) - DT50[[kinmodel]] = switch(kinmodel, - SFO = log(2)/coef(m)[["k"]], + TSS[[kinmodel]] = sum((m$model$parent - mean(m$model$parent))^2) + + DT50.o = switch(kinmodel, + SFO = log(2)/coef(m)[["k"]], FOMC = coef(m)[["beta"]] * (2^(1/coef(m)[["alpha"]]) - 1), - HS = optimize(f[[kinmodel]], c(0, max(kindata$t)), x=50)$minimum, - DFOP = optimize(f[[kinmodel]], c(0, max(kindata$t)), x=50)$minimum) - DT90[[kinmodel]] = switch(kinmodel, - SFO = log(10)/coef(m)[["k"]], - FOMC = coef(m)[["beta"]] * (10^(1/coef(m)[["alpha"]]) - 1), - HS = optimize(f[[kinmodel]], c(0, max(kindata$t)), x=90)$minimum, - DFOP = optimize(f[[kinmodel]], c(0, max(kindata$t)), x=90)$minimum) - err.min[[kinmodel]] <- kinerrmin(kinfits, kinmodel) + HS = optimize(f[[kinmodel]], c(0, DTmax), x=50)$minimum, + DFOP = optimize(f[[kinmodel]], c(0, DTmax), x=50)$minimum) + + + DT50[[kinmodel]] = ifelse(abs(DT50.o - DTmax) < 0.1, NA, DT50.o) + DT90.o = switch(kinmodel, + SFO = log(10)/coef(m)[["k"]], + FOMC = coef(m)[["beta"]] * (10^(1/coef(m)[["alpha"]]) - 1), + HS = optimize(f[[kinmodel]], c(0, DTmax), x=90)$minimum, + DFOP = optimize(f[[kinmodel]], c(0, DTmax), x=90)$minimum) + DT90[[kinmodel]] = ifelse(abs(DT90.o - DTmax) < 0.1, NA, DT90.o) + + # Chi2 error level as defined in FOCUS kinetics (see ?kinerrmin) + err.min[[kinmodel]] <- kinerrmin(kinfits, kinmodel) + + # Coefficient of determination calculated from residual sum of squares and totals sum of squares + # so this r2 is what is called model efficiency in FOCUS kinetics (2006), p. 99 + R2[[kinmodel]] = 1 - RSS[[kinmodel]]/TSS[[kinmodel]] + } } - stats <- data.frame(n.times = n.times, df = df, mean.means = kindata.means.mean, - RSS = RSS, err.min = err.min) + + stats <- data.frame(n.times = n.times, df = df, + mean.means = kindata.means.mean, + RSS = RSS, err.min = err.min, R2 = R2) results <- data.frame(DT50 = DT50, DT90 = DT90) list(parms = parms, stats = stats, results = results) } diff --git a/R/kinwrite.KinGUI.R b/R/kinwrite.KinGUI.R index bf94f49..2445685 100644 --- a/R/kinwrite.KinGUI.R +++ b/R/kinwrite.KinGUI.R @@ -1,3 +1,23 @@ +# $Id: kinwrite.KinGUI.R 59 2010-07-28 12:29:15Z jranke $ + +# Copyright (C) 2008-2010 Johannes Ranke +# Contact: mkin-devel@lists.berlios.de + +# This file is part of the R package kinfit + +# kinfit 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 + kinwrite.KinGUI <- function(kinobject, file, comment=NA) { sink(file) -- cgit v1.2.1