From 30cbb4092f6d2d3beff5800603374a0d009ad770 Mon Sep 17 00:00:00 2001 From: jranke Date: Tue, 11 May 2010 23:03:37 +0000 Subject: Initial upload of the upcoming multicompartment version of kinfit. Some functionality is still missing (chi2), some may never be implemented (FOMC model), but in general it is much more powerful than kinfit, owing to the powerful FME package. git-svn-id: svn+ssh://svn.r-forge.r-project.org/svnroot/kinfit/pkg/mkin@8 edb9625f-4e0d-4859-8d74-9fd3b1da38cb --- DESCRIPTION | 13 ++++ R/mkinfit.R | 66 ++++++++++++++++++ R/mkinmod.R | 76 ++++++++++++++++++++ data/FOCUS_2006_datasets.RData | Bin 0 -> 905 bytes inst/doc/header.tex | 31 +++++++++ inst/doc/mkin.Rnw | 155 +++++++++++++++++++++++++++++++++++++++++ inst/doc/mkin.pdf | Bin 0 -> 172301 bytes inst/doc/references.bib | 55 +++++++++++++++ inst/doc/run.bat | 5 ++ inst/unitTests/Makefile | 15 ++++ inst/unitTests/runit.mkinmod.R | 53 ++++++++++++++ man/FOCUS_2006_datasets.Rd | 35 ++++++++++ man/mkinfit.Rd | 85 ++++++++++++++++++++++ man/mkinmod.Rd | 44 ++++++++++++ tests/doRUnit.R | 58 +++++++++++++++ 15 files changed, 691 insertions(+) create mode 100644 DESCRIPTION create mode 100644 R/mkinfit.R create mode 100644 R/mkinmod.R create mode 100644 data/FOCUS_2006_datasets.RData create mode 100644 inst/doc/header.tex create mode 100644 inst/doc/mkin.Rnw create mode 100644 inst/doc/mkin.pdf create mode 100644 inst/doc/references.bib create mode 100644 inst/doc/run.bat create mode 100644 inst/unitTests/Makefile create mode 100644 inst/unitTests/runit.mkinmod.R create mode 100644 man/FOCUS_2006_datasets.Rd create mode 100644 man/mkinfit.Rd create mode 100644 man/mkinmod.Rd create mode 100644 tests/doRUnit.R diff --git a/DESCRIPTION b/DESCRIPTION new file mode 100644 index 0000000..96f651e --- /dev/null +++ b/DESCRIPTION @@ -0,0 +1,13 @@ +Package: mkin +Type: Package +Title: Routines for fitting kinetic models with one or more state variables to chemical degradation data +Version: 0.5-1 +Date: 2010-05-11 +Author: Johannes Ranke +Maintainer: Johannes Ranke +Description: Calculation routines based on the FOCUS Kinetics + Report (2006) +Depends: FME +License: GPL +LazyLoad: yes +LazyData: yes diff --git a/R/mkinfit.R b/R/mkinfit.R new file mode 100644 index 0000000..9651fd6 --- /dev/null +++ b/R/mkinfit.R @@ -0,0 +1,66 @@ +mkinfit <- function(mkinmod, observed, + parms.ini = rep(0.1, length(mkinmod$parms)), + state.ini = c(100, rep(0, length(mkinmod$diffs) - 1)), + fixed_parms = rep(FALSE, length(mkinmod$parms)), + fixed_initials = c(FALSE, rep(TRUE, length(mkinmod$diffs) - 1)), + plot = NULL, + err = NULL, weight = "none", scaleVar = FALSE, + ...) +{ + # Name the parameters if they are not named yet + if(is.null(names(parms.ini))) names(parms.ini) <- mkinmod$parms + # Create a function calculating the differentials specified by the model + mkindiff <- function(t, state, parms) { + diffs <- vector() + for (box in names(mkinmod$diffs)) + { + diffname <- paste("d", box, sep="_") + diffs[diffname] <- with(as.list(c(state, parms)), + eval(parse(text=mkinmod$diffs[[box]]))) + } + return(list(c(diffs))) + } + + # Name the inital parameter values if they are not named yet + if(is.null(names(state.ini))) names(state.ini) <- names(mkinmod$diffs) + + # TODO: Collect parameters to be optimised + parms.optim <- parms.ini[!fixed_parms] + parms.fixed <- parms.ini[fixed_parms] + + state.ini.optim <- state.ini[!fixed_initials] + state.ini.optim.boxnames <- names(state.ini.optim) + names(state.ini.optim) <- paste(names(state.ini.optim), "0", sep="_") + state.ini.fixed <- state.ini[fixed_initials] + + # Define the model cost function + cost <- function(P) + { + if(length(state.ini.optim) > 0) { + odeini <- c(P[1:length(state.ini.optim)], state.ini.fixed) + names(odeini) <- c(state.ini.optim.boxnames, names(state.ini.fixed)) + } else odeini <- state.ini.fixed + + odeparms <- c(P[(length(state.ini.optim) + 1):length(P)], parms.fixed) + # Solve the ODE + out <- ode( + y = odeini, + times = unique(observed$time), + func = mkindiff, + parms = odeparms) + + # Output transformation for models with ghost compartments like SFORB + out_transformed <- data.frame(time = out[,"time"]) + for (var in names(mkinmod$map)) { + if(length(mkinmod$map[[var]]) == 1) { + out_transformed[var] <- out[, var] + } else { + out_transformed[var] <- rowSums(out[, mkinmod$map[[var]]]) + } + } + + return(modCost(out_transformed, observed, y = "value", + err = err, weight = weight, scaleVar = scaleVar)) + } + modFit(cost, c(state.ini.optim, parms.optim), ...) +} diff --git a/R/mkinmod.R b/R/mkinmod.R new file mode 100644 index 0000000..7e64f87 --- /dev/null +++ b/R/mkinmod.R @@ -0,0 +1,76 @@ +mkinmod <- function(spec = list(parent = list(type = "SFO", to = NA, sink = TRUE))) +{ + if(!is.list(spec)) stop("spec must be a list of model specifications for each observed variable") + + # The returned model will be a list of two character vectors, containing + # differential equations and parameter names + parms <- vector() + diffs <- vector() + map <- list() + + # Establish list of differential equations + for (varname in names(spec)) + { + if(!spec[[varname]]$type %in% c("SFO", "SFORB")) stop("Available types are SFO and SFORB only") + new_parms <- vector() + + # New (sub)compartments (boxes) needed for the model type + new_boxes <- switch(spec[[varname]]$type, + SFO = varname, + SFORB = paste(varname, c("free", "bound"), sep="_") + ) + map[[varname]] <- new_boxes + + # Start a new differential equation for each new box + new_diffs <- paste("d_", new_boxes, " =", sep="") + + # Construct terms for transfer to sink and add if appropriate + if(spec[[varname]]$sink) { + k_compound_sink <- paste("k", new_boxes[[1]], "sink", sep="_") + sink_term <- paste("-", k_compound_sink, "*", new_boxes[[1]]) + # Only add sink term to first (or only) box for SFO and SFORB models + if(spec[[varname]]$type %in% c("SFO", "SFORB")) { + new_diffs[[1]] <- paste(new_diffs[[1]], sink_term) + } + new_parms <- k_compound_sink + } + + # Add reversible binding if appropriate + if(spec[[varname]]$type == "SFORB") { + k_free_bound <- paste("k", varname, "free", "bound", sep="_") + k_bound_free <- paste("k", varname, "bound", "free", sep="_") + reversible_binding_terms <- c( + paste("-", k_free_bound, "*", new_boxes[[1]], "+", k_bound_free, "*", new_boxes[[2]]), + paste("+", k_free_bound, "*", new_boxes[[1]], "-", k_bound_free, "*", new_boxes[[2]])) + new_diffs <- paste(new_diffs, reversible_binding_terms) + new_parms <- c(new_parms, k_free_bound, k_bound_free) + } + + # Add observed variable to model + parms <- c(parms, new_parms) + names(new_diffs) <- new_boxes + diffs <- c(diffs, new_diffs) + } + + # Transfer between compartments + for (varname in names(spec)) { + to <- spec[[varname]]$to + if(!is.na(to)) { + origin_box <- switch(spec[[varname]]$type, + SFO = varname, + SFORB = paste(varname, "free", sep="_")) + for (target in to) { + target_box <- switch(spec[[target]]$type, + SFO = target, + SFORB = paste(target, "free", sep="_")) + } + k_from_to <- paste("k", origin_box, target_box, sep="_") + diffs[[origin_box]] <- paste(diffs[[origin_box]], "-", k_from_to, "*", origin_box) + diffs[[target_box]] <- paste(diffs[[target_box]], "+", k_from_to, "*", origin_box) + parms <- c(parms, k_from_to) + } + } + model <- list(diffs = diffs, parms = parms, map = map) + class(model) <- "mkinmod" + return(model) +} diff --git a/data/FOCUS_2006_datasets.RData b/data/FOCUS_2006_datasets.RData new file mode 100644 index 0000000..eb9e526 Binary files /dev/null and b/data/FOCUS_2006_datasets.RData differ diff --git a/inst/doc/header.tex b/inst/doc/header.tex new file mode 100644 index 0000000..9d6ec49 --- /dev/null +++ b/inst/doc/header.tex @@ -0,0 +1,31 @@ +\usepackage{booktabs} +\usepackage{amsfonts} +\usepackage{latexsym} +\usepackage{amsmath} +\usepackage{amssymb} +\usepackage{graphicx} +\usepackage{parskip} +\usepackage[round]{natbib} +\usepackage{amstext} +\usepackage{hyperref} +\usepackage[latin1]{inputenc} + +\newcommand{\Rpackage}[1]{{\normalfont\fontseries{b}\selectfont #1}} +\newcommand{\Robject}[1]{\texttt{#1}} +\newcommand{\Rclass}[1]{\textit{#1}} +\newcommand{\Rcmd}[1]{\texttt{#1}} + +\newcommand{\RR}{\textsf{R}} + +\RequirePackage[T1]{fontenc} +\RequirePackage{graphicx,ae,fancyvrb} +\IfFileExists{upquote.sty}{\RequirePackage{upquote}}{} +\usepackage{relsize} + +\DefineVerbatimEnvironment{Sinput}{Verbatim}{baselinestretch=1.05} +\DefineVerbatimEnvironment{Soutput}{Verbatim}{fontfamily=courier, + baselinestretch=1.05, + fontshape=it, + fontsize=\relsize{-1}} +\DefineVerbatimEnvironment{Scode}{Verbatim}{} +\newenvironment{Schunk}{}{} diff --git a/inst/doc/mkin.Rnw b/inst/doc/mkin.Rnw new file mode 100644 index 0000000..84ede0f --- /dev/null +++ b/inst/doc/mkin.Rnw @@ -0,0 +1,155 @@ +%%\VignetteIndexEntry{Routines for fitting kinetic models with one or more state variables to chemical degradation data} +%%VignetteDepends{FME} +%%\usepackage{Sweave} +\documentclass[12pt,a4paper]{article} +\usepackage{a4wide} +%%\usepackage[lists,heads]{endfloat} +\input{header} +\hypersetup{ + pdftitle = {mkin - Routines for fitting kinetic models with one or more state variables to chemical degradation data}, + pdfsubject = {Manuscript}, + pdfauthor = {Johannes Ranke}, + colorlinks = {true}, + linkcolor = {blue}, + citecolor = {blue}, + urlcolor = {red}, + hyperindex = {true}, + linktocpage = {true}, +} +\SweaveOpts{engine=R, eps=FALSE, keep.source = TRUE} +<>= +options(prompt = "R> ") +options(SweaveHooks = list( + cex = function() par(cex.lab = 1.3, cex.axis = 1.3))) +@ +\begin{document} +\title{mkin -\\ +Routines for fitting kinetic models with one or more state variables to chemical degradation data} +\author{\textbf{Johannes Ranke} \\ +%EndAName +Product Safety \\ +Harlan Laboratories Ltd. \\ +Zelgliweg 1, CH--4452 Itingen, Switzerland} +\maketitle + +\begin{abstract} +In the regulatory evaluation of chemical substances like plant protection +products (pesticides), biocides and other chemicals, degradation data play an +important role. For the evaluation of pesticide degradation experiments, +detailed guidance has been developed, based on nonlinear optimisation. +The \RR{} add-on package \Rpackage{mkin} implements fitting some of the models +recommended in this guidance from within R and calculates some statistical +measures for data series within one or more compartments, for parent and +metabolites. +\end{abstract} + + +\thispagestyle{empty} \setcounter{page}{0} + +\clearpage + +\tableofcontents + +\textbf{Key words}: Kinetics, FOCUS, nonlinear optimisation + +\section{Introduction} +\label{intro} + +Many approaches are possible regarding the evaluation of chemical degradation +data. The \Rpackage{kinfit} package \citep{pkg:kinfit} in \RR{} +\citep{rcore2010} implements the approach recommended in the kinetics report +provided by the FOrum for Co-ordination of pesticide fate models and their +USe \citep{FOCUS2006} for simple data series for one parent compound in one +compartment. + +The \Rpackage{mkin} package \citep{pkg:mkin} extends this approach to data series +with metabolites and more than one compartment and includes the possibility +for back reactions. + +\section{Example} +\label{exam} + +In the following, requirements for data formatting are explained. Then the +procedure for fitting the four kinetic models recommended by the FOCUS group +to an example dataset for parent only given in the FOCUS kinetics report is +illustrated. The explanations are kept rather verbose in order to lower the +barrier for \RR{} newcomers. + +\subsection{Data format} + +The following listing shows example dataset C from the FOCUS kinetics +report as distributed with the \Rpackage{kinfit} package + +<>= +library("mkin") +FOCUS_2006_C +@ + +Note that the data needs to be in the format of a data frame containing a +variable \Robject{name} specifying the observed variable, indicating the +compound name and, if applicable, the compartment, a variable \Robject{time} +containing sampling times, and a numeric variable \Robject{value} specifying +the observed value of the variable. If a further variable \Robject{error} +is present, this will be used to give different weights to the data points +(the higher the error, the lower the weight, see the help page of the +\Robject{modCost} function of the \Rpackage{FME} package \citep{soetaert10}). +Replicate measurements are not recorded in extra columns but simply appended, +leading to multiple occurrences of the sampling times \Robject{time}. + +Small to medium size dataset can be conveniently entered directly as \RR{} code +as shown in the following listing + +<>= +example_data <- data.frame( + time = c(0, 1, 3, 7, 14, 28, 63, 91, 119), + parent = c(85.1, 57.9, 29.9, 14.6, 9.7, 6.6, 4, 3.9, 0.6) +) +@ + +\subsection{Model definition} + +The next task is to define the model to be fitted to the data. In order to +facilitate this task, a convenience function \Robject{mkinmod} is available. + +<>= +SFO <- mkinmod(spec = list(parent = list(type = "SFO", to = NA, sink = TRUE))) +SFORB <- mkinmod(spec = list(parent = list(type = "SFORB", to = NA, sink = TRUE))) +SFO_SFO <- mkinmod(spec = list( + parent = list(type = "SFO", to = "m1", sink = TRUE), + m1 = list(type = "SFO", to = NA, sink = TRUE))) +SFORB_SFO <- mkinmod(spec = list( + parent = list(type = "SFORB", to = "m1", sink = TRUE), + m1 = list(type = "SFO", to = NA, sink = TRUE))) +@ + +\subsection{Fitting the model} + +Then the model parameters should be fitted to the data. The function +\Robject{mkinfit} internally creates a cost function using \Robject{modCost} +from the \Rpackage{FME} package and the produces a fit using \Robject{modFit} +from the same package. + +<>= +# Do not show significance stars as they interfere with vignette generation +options(show.signif.stars = FALSE) +SFO.fit <- mkinfit(SFO, FOCUS_2006_C) +summary(SFO.fit) +SFORB.fit <- mkinfit(SFORB, FOCUS_2006_C) +summary(SFORB.fit) +SFO_SFO.fit <- mkinfit(SFO_SFO, FOCUS_2006_D) +summary(SFO_SFO.fit) +SFO_SFO.fit.2 <- mkinfit(SFO_SFO, FOCUS_2006_D, + fixed_initials = c(FALSE, FALSE), fixed_parms = c(FALSE, TRUE, FALSE)) +summary(SFO_SFO.fit.2) +SFO_SFO.fit.3 <- mkinfit(SFO_SFO, FOCUS_2006_D, + fixed_initials = c(FALSE, FALSE), fixed_parms = c(FALSE, TRUE, FALSE), lower = -0.0000001) +summary(SFO_SFO.fit.3) +SFORB_SFO.fit <- mkinfit(SFORB_SFO, FOCUS_2006_D) +summary(SFORB_SFO.fit) +@ + +\bibliographystyle{plainnat} +\bibliography{references} + +\end{document} +% vim: set foldmethod=syntax: diff --git a/inst/doc/mkin.pdf b/inst/doc/mkin.pdf new file mode 100644 index 0000000..95c8a69 Binary files /dev/null and b/inst/doc/mkin.pdf differ diff --git a/inst/doc/references.bib b/inst/doc/references.bib new file mode 100644 index 0000000..b158eeb --- /dev/null +++ b/inst/doc/references.bib @@ -0,0 +1,55 @@ +@Manual{pkg:kinfit, + title = {kinfit: {R}outines for fitting simple kinetic models to chemical degradation data}, + author = {Johannes Ranke}, + year = {2010}, + url = {http://CRAN.R-project.org} +} + +@Manual{pkg:mkin, + title = {mkin: {R}outines for fitting kinetic models with one or more state variables to chemical degradation data}, + author = {Johannes Ranke}, + year = {2010}, + url = {http://CRAN.R-project.org} +} + +@Article{soetaert10, + title = {Inverse Modelling, Sensitivity and Monte Carlo Analysis in {R} Using Package {FME}}, + author = {Karline Soetaert and Thomas Petzoldt}, + journal = {Journal of Statistical Software}, + year = {2010}, + volume = {33}, + number = {3}, + pages = {1--28}, + url = {http://www.jstatsoft.org/v33/i03/} +} + +@Manual{ rcore2010, + title = {\textsf{R}: A Language and Environment for Statistical + Computing}, + author = {{R Development Core Team}}, + organization = {R Foundation for Statistical Computing}, + address = {Vienna, Austria}, + year = 2010, + note = {{ISBN} 3-900051-07-0}, + url = {http://www.R-project.org} +} + +@Manual{ FOCUS2006, + title = {Guidance Document on Estimating Persistence and + Degradation Kinetics from Environmental Fate Studies on Pesticides in EU + Registration. Report of the FOCUS Work Group on Degradation Kinetics}, + note = {EC Document Reference Sanco/10058/2005 version 2.0}, + author = {{FOCUS Work Group on Degradation Kinetics}}, + year = {2006}, + url = {http://focus.jrc.ec.europa.eu/dk} +} + +@Inproceedings{ schaefer2007, + title = {{KinGUI}: a new kinetic software tool for evaluations according to FOCUS degradation kinetics}, + author = {D. Sch\"{a}fer and M. Mikolasch and P. Rainbird and B. Harvey}, + booktitle = {Proceedings of the XIII Symposium Pesticide Chemistry}, + editor = {Del Re A. A. M. and Capri E. and Fragoulis G. and Trevisan M.}, + year = {2007}, + address = {Piacenza}, + pages = {916--923} +} diff --git a/inst/doc/run.bat b/inst/doc/run.bat new file mode 100644 index 0000000..c28c666 --- /dev/null +++ b/inst/doc/run.bat @@ -0,0 +1,5 @@ +R.exe -e "Sweave('mkin.Rnw', stylepath=FALSE)" +pdflatex.exe mkin +bibtex.exe mkin +pdflatex.exe mkin +pdflatex.exe mkin diff --git a/inst/unitTests/Makefile b/inst/unitTests/Makefile new file mode 100644 index 0000000..8d13225 --- /dev/null +++ b/inst/unitTests/Makefile @@ -0,0 +1,15 @@ +TOP=../.. +PKG=${shell cd ${TOP};pwd} +SUITE=doRUnit.R +R=R + +all: inst test + +inst: # Install package + cd ${TOP}/..;\ + ${R} CMD INSTALL ${PKG} + +test: # Run unit tests + export RCMDCHECK=FALSE;\ + cd ${TOP}/tests;\ + ${R} --vanilla --slave < ${SUITE} diff --git a/inst/unitTests/runit.mkinmod.R b/inst/unitTests/runit.mkinmod.R new file mode 100644 index 0000000..b6ca6b8 --- /dev/null +++ b/inst/unitTests/runit.mkinmod.R @@ -0,0 +1,53 @@ +test.mkinmod.SFO <- function() +{ + SFO.diffs <- c( + parent = "d_parent = - k_parent_sink * parent" + ) + SFO.parms <- c("k_parent_sink") + SFO.map <- list(parent = "parent") + SFO <- list(diffs = SFO.diffs, parms = SFO.parms, map = SFO.map) + class(SFO) <- "mkinmod" + SFO.mkinmod <- mkinmod(spec = list( + parent = list(type = "SFO", to = NA, sink=TRUE)) + ) + checkIdentical(SFO, SFO.mkinmod) +} + +test.mkinmod.SFORB <- function() +{ + SFORB.diffs <- c( + parent_free = paste( + "d_parent_free = - k_parent_free_sink * parent_free", + "- k_parent_free_bound * parent_free", + "+ k_parent_bound_free * parent_bound"), + parent_bound = paste( + "d_parent_bound =", + "+ k_parent_free_bound * parent_free", + "- k_parent_bound_free * parent_bound") + ) + SFORB.parms <- c("k_parent_free_sink", "k_parent_free_bound", "k_parent_bound_free") + SFORB.map <- list(parent = c("parent_free", "parent_bound")) + SFORB <- list(diffs = SFORB.diffs, parms = SFORB.parms, map = SFORB.map) + class(SFORB) <- "mkinmod" + SFORB.mkinmod <- mkinmod(spec = list( + parent = list(type = "SFORB", to = NA, sink=TRUE)) + ) + checkIdentical(SFORB, SFORB.mkinmod) +} + +test.mkinmod.SFO_SFO <- function() +{ + SFO_SFO.diffs <- c( + parent = "d_parent = - k_parent_sink * parent - k_parent_m1 * parent", + m1 = "d_m1 = - k_m1_sink * m1 + k_parent_m1 * parent" + ) + SFO_SFO.parms <- c("k_parent_sink", "k_m1_sink", "k_parent_m1") + SFO_SFO.map <- list(parent = "parent", m1 = "m1") + SFO_SFO <- list(diffs = SFO_SFO.diffs, parms = SFO_SFO.parms, map = SFO_SFO.map) + class(SFO_SFO) <- "mkinmod" + SFO_SFO.mkinmod <- mkinmod(spec = list( + parent = list(type = "SFO", to = "m1", sink=TRUE), + m1 = list(type = "SFO", to = NA, sink=TRUE)) + ) + checkIdentical(SFO_SFO, SFO_SFO.mkinmod) +} diff --git a/man/FOCUS_2006_datasets.Rd b/man/FOCUS_2006_datasets.Rd new file mode 100644 index 0000000..04565f7 --- /dev/null +++ b/man/FOCUS_2006_datasets.Rd @@ -0,0 +1,35 @@ +\name{FOCUS_2006_datasets} +\Rdversion{1.1} +\alias{FOCUS_2006_A} +\alias{FOCUS_2006_B} +\alias{FOCUS_2006_C} +\alias{FOCUS_2006_D} +\alias{FOCUS_2006_E} +\alias{FOCUS_2006_F} +\docType{data} +\title{ +Datasets A to F from the FOCUS Kinetics report from 2006 +} +\description{ +Data taken from an FOCUS (2006), p. 258. +} +\usage{FOCUS_2006_datasets} +\format{ + 6 datasets with observations on the following variables. + \describe{ + \item{\code{name}}{a factor containing the name of the observed variable} + \item{\code{time}}{a numeric vector containing time points} + \item{\code{value}}{a numeric vector containing concentrations in percent of applied radioactivity} + } +} +\source{ + FOCUS (2006) \dQuote{Guidance Document on Estimating Persistence and + Degradation Kinetics from Environmental Fate Studies on Pesticides in EU + Registration} Report of the FOCUS Work Group on Degradation Kinetics, + EC Document Reference Sanco/10058/2005 version 2.0, 434 pp, + \url{http://focus.jrc.ec.europa.eu/dk} +} +\examples{ +FOCUS_2006_C +} +\keyword{datasets} diff --git a/man/mkinfit.Rd b/man/mkinfit.Rd new file mode 100644 index 0000000..7ddc10c --- /dev/null +++ b/man/mkinfit.Rd @@ -0,0 +1,85 @@ +\name{mkinfit} +\alias{mkinfit} +\title{ + Fit a kinetic model to data with one or more state variables. +} +\description{ + This function uses the Flexible Modelling Environment package + \code{\link{FME}} to create a function calculating the model cost, which is + then minimised, using the specified initial or fixed parameters and starting + values. +} +\usage{ +mkinfit(mkinmod, observed, parms.ini = rep(0.1, length(mkinmod$parms)), state.ini = c(100, rep(0, length(mkinmod$diffs) - 1)), fixed_parms = rep(FALSE, length(mkinmod$parms)), fixed_initials = c(FALSE, rep(TRUE, length(mkinmod$diffs) - 1)), plot = NULL, err = NULL, weight = "none", scaleVar = FALSE, ...) +} +\arguments{ + \item{mkinmod}{ + A list of class \code{\link{mkinmod}}, containing the kinetic model to be fitted to the data. + } + \item{observed}{ + The observed data. It has to be in the long format as described in + \code{\link{modFit}}, i.e. the first column must contain the name of the + observed variable for each data point. The second column must contain the + times of observation, named "time". The third column must be named "value" + and contain the observed values. Optionally, a further column can contain + weights for each data point. If it is not named "err", its name must be + passed as a further argument named \code{err} which is then passed on to + \code{\link{modFit}}. + } + \item{parms.ini}{ + A named vector if initial values for the parameters, including both parameters to + be optimised and potentially also fixed parameters as indicated by \code{fixed_parms}. +} + \item{state.ini}{ + A named vector of initial values for the state variables of the model. In case the + observed variables are represented by more than one model variable, the names will + differ from the names of the observed variables (see \code{map} component of + \code{\link{mkinmod}}). The default is to set the initial value of the first model + variable to 100 and all others to 0. +} + \item{fixed_parms}{ + A vector of booleans specifying which parameters are not to be optimised. The default + is to include all model parameters in the optimisation. +} + \item{fixed_initials}{ + A vector of booleans specifying which initial values to include in the optimisation. + The default is to optimise the initial value of the first model variable and to + keep all other initial values fixed. +} + \item{plot}{ + Should the observed values and the numerical solutions be plotted at each stage + of the optimisation? +} + \item{err }{either \code{NULL}, or the name of the column with the + \emph{error} estimates, used to weigh the residuals (see details of + \code{\link{modCost}}); if \code{NULL}, then the residuals are not weighed. +} + \item{weight}{only if \code{err}=\code{NULL}: how to weigh the + residuals, one of "none", "std", "mean", see details of \code{\link{modCost}}. +} + \item{scaleVar}{ + Will be passed to \code{\link{modCost}}. Default is not to scale Variables according + to the number of observations. +} + \item{\dots}{ + Further arguments that will be passed to \code{\link{modFit}}. +} +} +\value{ + A list of class \code{\link{modFit}}. Thus, at present, a summary can be obtained + by \code{\link{summary.modFit}}. +} +\author{ + Johannes Ranke +} +\examples{ +# One parent compound, one metabolite, both single first order. +SFO_SFO <- mkinmod(spec = list( + parent = list(type = "SFO", to = "m1", sink = TRUE), + m1 = list(type = "SFO", to = NA, sink = TRUE))) +# Fit the model to the FOCUS example dataset D using defaults +fit <- mkinfit(SFO_SFO, FOCUS_2006_D) +summary(fit) +} +\keyword{ models } +\keyword{ optimize } diff --git a/man/mkinmod.Rd b/man/mkinmod.Rd new file mode 100644 index 0000000..4569b31 --- /dev/null +++ b/man/mkinmod.Rd @@ -0,0 +1,44 @@ +\name{mkinmod} +\alias{mkinmod} +\title{ + Function to set up a kinetic model with one or more state variables. +} +\description{ + The function takes a specification, consisting of a list of the observed variables + in the data. Each observed variable is again represented by a list, specifying the + kinetic model type and reaction or transfer to other observed compartments. +} +\usage{ +mkinmod(spec = list(parent = list(type = "SFO", to = NA, sink = TRUE))) +} +\arguments{ + \item{spec}{ + A list of observed variables to be modelled. Each observed variable has to be + represented by a list with the following entries: + \code{type}{ The type of kinetics to use for the variable. Currently, only + single first order kinetics "SFO" or single first order with reversible binding + "SFORB" are implemented. } + \code{to}{ A vector of names of variables to which a transfer is to be assumed + in the model. } + \code{sink}{ Boolean, specifying if transformation to unspecified compounds (sink) + is to be assumed in the model. } + } +} +\value{ + A list of class \code{mkinmod} for use with \code{\link{mkinfit}}, containing + \item{diffs}{ A vector of string representations of differential equations, + one for each modelling variable. } + \item{parms}{ A vector of parameter names occurring in the differential equations. } + \item{map}{ A list containing named character vectors for each observed variable, specifying + the modelling variables by which it is represented. } +} +\author{ + Johannes Ranke +} +\examples{ +# One parent compound, one metabolite, both single first order. +SFO_SFO <- mkinmod(spec = list( + parent = list(type = "SFO", to = "m1", sink = TRUE), + m1 = list(type = "SFO", to = NA, sink = TRUE))) +} +\keyword{ models } diff --git a/tests/doRUnit.R b/tests/doRUnit.R new file mode 100644 index 0000000..e2ad4cd --- /dev/null +++ b/tests/doRUnit.R @@ -0,0 +1,58 @@ +if(require("RUnit", quietly=TRUE)) { + + ## --- Setup --- + + pkg <- "mkin" # <-- Change to package name! + if(Sys.getenv("RCMDCHECK") == "FALSE") { + ## Path to unit tests for standalone running under Makefile (not R CMD check) + ## PKG/tests/../inst/unitTests + path <- file.path(getwd(), "..", "inst", "unitTests") + } else { + ## Path to unit tests for R CMD check + ## PKG.Rcheck/tests/../PKG/unitTests + path <- system.file(package=pkg, "unitTests") + } + cat("\nRunning unit tests\n") + print(list(pkg=pkg, getwd=getwd(), pathToUnitTests=path)) + + library(package=pkg, character.only=TRUE) + + ## If desired, load the name space to allow testing of private functions + ## if (is.element(pkg, loadedNamespaces())) + ## attach(loadNamespace(pkg), name=paste("namespace", pkg, sep=":"), pos=3) + ## + ## or simply call PKG:::myPrivateFunction() in tests + + ## --- Testing --- + + ## Define tests + testSuite <- defineTestSuite(name=paste(pkg, " Unit Tests"), + dirs=path) + ## Run + tests <- runTestSuite(testSuite) + + ## Default report name + pathReport <- file.path(path, "report") + + ## Report to stdout and text files + cat("------------------- UNIT TEST SUMMARY ---------------------\n\n") + printTextProtocol(tests, showDetails=FALSE) + printTextProtocol(tests, showDetails=FALSE, + fileName=paste(pathReport, "Summary.txt", sep="")) + printTextProtocol(tests, showDetails=TRUE, + fileName=paste(pathReport, ".txt", sep="")) + + ## Report to HTML file + printHTMLProtocol(tests, fileName=paste(pathReport, ".html", sep="")) + + ## Return stop() to cause R CMD check stop in case of + ## - failures i.e. FALSE to unit tests or + ## - errors i.e. R errors + tmp <- getErrors(tests) + if(tmp$nFail > 0 | tmp$nErr > 0) { + stop(paste("\n\nunit testing failed (#test failures: ", tmp$nFail, + ", #R errors: ", tmp$nErr, ")\n\n", sep="")) + } +} else { + warning("cannot run unit tests -- package RUnit is not available") +} -- cgit v1.2.1