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 --- ChangeLog | 12 ++++++++++-- DESCRIPTION | 2 +- R/mkinerrmin.R | 15 +++++++++------ 3 files changed, 20 insertions(+), 9 deletions(-) diff --git a/ChangeLog b/ChangeLog index a57894d3..422deb2c 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,3 +1,9 @@ +2013-10-21 Johannes Ranke for mkin (0.9-22) + + * Get rid of the optimisation step in mkinerrmin - this was unnecessary + Thanks to KinGUII for the inspiration - actually this is equation 6-2 + in FOCUS kinetics p. 91 that I had overlooked originally + 2013-10-17 Johannes Ranke for mkin (0.9-22) * Fix plot.mkinfit as it passed graphical arguments like main to the solver @@ -16,11 +22,13 @@ * Correct the output of the data in the case of manual weighting * Implement IRLS assuming different variances for observed variables - 2013-10-09 Johannes Ranke for version (0.9-22) * Do not use 0 values at time zero for chi2 error level calculations. - This is the way it is done in KinGUII and it makes sense + This is the way it is done in KinGUII and it makes sense. It does + impact the chi2 error levels in the output. Generally they seem to be + lower for metabolites now, presumably because the mean of the observed + values is higher Changes performed in earlier versions are documented in the subversion log files on R-Forge http://www.r-forge.r-project.org/scm/?group_id=615 diff --git a/DESCRIPTION b/DESCRIPTION index f7a293d3..9e003570 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -3,7 +3,7 @@ Type: Package Title: Routines for fitting kinetic models with one or more state variables to chemical degradation data Version: 0.9-22 -Date: 2013-10-10 +Date: 2013-10-21 Author: Johannes Ranke, with contributions from Katrin Lindenberger, René Lehmann Maintainer: Johannes Ranke Description: Calculation routines based on the FOCUS Kinetics Report (2006). diff --git a/R/mkinerrmin.R b/R/mkinerrmin.R index 7d9d4626..074cc568 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