15 files changed, 691 insertions, 0 deletions
diff --git a/DESCRIPTION b/DESCRIPTION
new file mode 100644
index 00000000..96f651eb
--- /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 <jranke@harlan.com>
+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 00000000..9651fd66
--- /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 00000000..7e64f879
--- /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 00000000..eb9e5266
--- /dev/null
+++ b/data/FOCUS_2006_datasets.RData
diff --git a/inst/doc/header.tex b/inst/doc/header.tex
new file mode 100644
index 00000000..9d6ec49b
--- /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 00000000..84ede0f5
--- /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}
+<<setup, echo = FALSE, results = hide>>=
+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
+
+<<FOCUS_2006_C_data, echo=TRUE, eval=TRUE>>=
+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
+
+<<data_format, echo=TRUE>>=
+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.
+
+<<model_definition, echo=TRUE>>=
+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.
+
+<<model_fitting, echo=TRUE>>=
+# 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 00000000..95c8a690
--- /dev/null
+++ b/inst/doc/mkin.pdf
diff --git a/inst/doc/references.bib b/inst/doc/references.bib
new file mode 100644
index 00000000..b158eebf
--- /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 00000000..c28c6663
--- /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 00000000..8d132253
--- /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 00000000..b6ca6b8d
--- /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 00000000..04565f74
--- /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 00000000..7ddc10c1
--- /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 <jranke@{harlan.com,uni-bremen.de}>
+}
+\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 00000000..4569b317
--- /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 <jranke@{harlan.com,uni-bremen.de}>
+}
+\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 00000000..e2ad4cd2
--- /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")
+}