From 9b5f45d5f45722875adada2afc554cec8308a761 Mon Sep 17 00:00:00 2001 From: jranke Date: Mon, 21 Oct 2013 09:46:27 +0000 Subject: - Simplified implementation of mkinerrmin (see ChangeLog) - Improve ChangeLog and add comment to source code - Simple test of the new mkinerrmin implementation in test.R - Added some summaries for documenting the changes in mkinerrmin git-svn-id: svn+ssh://svn.r-forge.r-project.org/svnroot/kinfit/pkg/mkin@120 edb9625f-4e0d-4859-8d74-9fd3b1da38cb --- R/mkinerrmin.R | 15 +++++++++------ 1 file changed, 9 insertions(+), 6 deletions(-) (limited to 'R') diff --git a/R/mkinerrmin.R b/R/mkinerrmin.R index 7d9d462..074cc56 100644 --- a/R/mkinerrmin.R +++ b/R/mkinerrmin.R @@ -22,16 +22,14 @@ if(getRversion() >= '2.15.1') utils::globalVariables(c("name", "value_mean")) mkinerrmin <- function(fit, alpha = 0.05) { parms.optim <- fit$par + kinerrmin <- function(errdata, n.parms) { means.mean <- mean(errdata$value_mean, na.rm = TRUE) df = length(errdata$value_mean) - n.parms - f <- function(err) - { - (sum((errdata$value_mean - errdata$value_pred)^2/((err * means.mean)^2)) - - qchisq(1 - alpha,df))^2 - } - err.min <- optimize(f, c(0.01,0.9))$minimum + err.min <- sqrt((1 / qchisq(1 - alpha, df)) * + sum((errdata$value_mean - errdata$value_pred)^2)/(means.mean^2)) + return(list(err.min = err.min, n.optim = n.parms, df = df)) } @@ -39,7 +37,12 @@ mkinerrmin <- function(fit, alpha = 0.05) errdata <- merge(means, fit$predicted, by = c("time", "name"), suffixes = c("_mean", "_pred")) errdata <- errdata[order(errdata$time, errdata$name), ] + + # Any value that is set to exactly zero is not really an observed value + # Remove those at time 0 - those are caused by the FOCUS recommendation + # to set metabolites occurring at time 0 to 0 errdata <- subset(errdata, !(time == 0 & value_mean == 0)) + n.optim.overall <- length(parms.optim) errmin.overall <- kinerrmin(errdata, n.optim.overall) -- cgit v1.2.1