diff options
Diffstat (limited to 'R')
-rw-r--r-- | R/D24_2014.R | 28 | ||||
-rw-r--r-- | R/f_time_norm_focus.R | 77 | ||||
-rw-r--r-- | R/focus_soil_moisture.R | 11 | ||||
-rw-r--r-- | R/mkin_wide_to_long.R | 2 | ||||
-rw-r--r-- | R/mkinds.R | 117 | ||||
-rw-r--r-- | R/mkinerrmin.R | 2 | ||||
-rw-r--r-- | R/mkinerrplot.R | 2 | ||||
-rw-r--r-- | R/mkinfit.R | 2 | ||||
-rw-r--r-- | R/mkinresplot.R | 2 | ||||
-rw-r--r-- | R/plot.mixed.mmkin.R | 2 | ||||
-rw-r--r-- | R/plot.mkinfit.R | 2 | ||||
-rw-r--r-- | R/saem.R | 2 |
12 files changed, 233 insertions, 16 deletions
diff --git a/R/D24_2014.R b/R/D24_2014.R new file mode 100644 index 00000000..81316833 --- /dev/null +++ b/R/D24_2014.R @@ -0,0 +1,28 @@ +#' Aerobic soil degradation data on 2,4-D from the EU assessment in 2014 +#' +#' The five datasets were extracted from the active substance evaluation dossier +#' published by EFSA. Kinetic evaluations shown for these datasets are intended +#' to illustrate and advance kinetic modelling. The fact that these data and +#' some results are shown here does not imply a license to use them in the +#' context of pesticide registrations, as the use of the data may be +#' constrained by data protection regulations. +#' +#' Metabolite residues at early sampling times reported as 0.0 were set to NA. +#' +#' The R code used to create this data object is installed with this package +#' in the 'dataset_generation' directory. In the code, page numbers are given for +#' specific pieces of information in the comments. +#' +#' @format An [mkindsg] object grouping five datasets +#' @source Hellenic Ministry of Rural Development and Agriculture (2014) +#' Final addendum to the Renewal Assessment Report - public version - 2,4-D +#' Volume 3 Annex B.8 Fate and behaviour in the environment p. 638, 640, +#' 644-646. +#' \url{http://registerofquestions.efsa.europa.eu/roqFrontend/outputLoader?output=ON-3812} +#' @examples +#' print(D24_2014) +#' print(D24_2014$ds[[1]], data = TRUE) +#' m1 = mkinmod(D24 = list(type = "SFO", to = "phenol"), +#' phenol = list(type = "SFO", to = "anisole"), +#' anisole = list(type = "SFO")) +"D24_2014" diff --git a/R/f_time_norm_focus.R b/R/f_time_norm_focus.R new file mode 100644 index 00000000..be5cf583 --- /dev/null +++ b/R/f_time_norm_focus.R @@ -0,0 +1,77 @@ +utils::globalVariables("D24_2014") + +#' Normalisation factors for aerobic soil degradation according to FOCUS guidance +#' +#' Time step normalisation factors for aerobic soil degradation as described +#' in Appendix 8 to the FOCUS kinetics guidance (FOCUS 2014, p. 369). +#' +#' @param object An object containing information used for the calculations +#' @param temperature Numeric vector of temperatures in °C +#' @param moisture Numeric vector of moisture contents in \\% w/w +#' @param field_moisture Numeric vector of moisture contents at field capacity +#' (pF2) in \\% w/w +#' @param Q10 The Q10 value used for temperature normalisation +#' @param walker The Walker exponent used for moisture normalisation +#' @param f_na The factor to use for NA values. If set to NA, only factors +#' for complete cases will be returned. +#' @param \dots Currently not used +#' @references +#' 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://esdac.jrc.ec.europa.eu/projects/degradation-kinetics} +#' FOCUS (2014) \dQuote{Generic guidance for Estimating Persistence +#' and Degradation Kinetics from Environmental Fate Studies on Pesticides in +#' EU Registration} Report of the FOCUS Work Group on Degradation Kinetics, +#' Version 1.1, 18 December 2014 +#' \url{http://esdac.jrc.ec.europa.eu/projects/degradation-kinetics} +#' @seealso [focus_soil_moisture] +#' @examples +#' f_time_norm_focus(25, 20, 25) # 1.37, compare p. 184 +#' +#' D24_2014$meta +#' # No moisture normalisation in the first dataset, so we use f_na = 1 to get +#' # Temperature only normalisation as in the EU evaluation +#' f_time_norm_focus(D24_2014, f_na = 1) +#' # Moisture normalisation for the other four soils is one, as soil moisture +#' # is higher than the approximate field capacity derived from the USDA soil +#' # type +#' @export +f_time_norm_focus <- function(object, ...) { + UseMethod("f_time_norm_focus") +} + +#' @rdname f_time_norm_focus +#' @export +f_time_norm_focus.numeric <- function(object, + moisture = NA, field_moisture = NA, + temperature = object, + Q10 = 2.58, walker = 0.7, f_na = NA, ...) +{ + f_temp <- ifelse(is.na(temperature), + f_na, + ifelse(temperature <= 0, + 0, + Q10^((temperature - 20)/10))) + f_moist <- ifelse(is.na(moisture), + f_na, + ifelse(moisture >= field_moisture, + 1, + (moisture / field_moisture)^walker)) + f_time_norm <- f_temp * f_moist + f_time_norm +} + +#' @rdname f_time_norm_focus +#' @export +f_time_norm_focus.mkindsg <- function(object, Q10 = 2.58, walker = 0.7, f_na = NA, ...) { + meta <- object$meta + field_moisture <- focus_soil_moisture[meta$usda_soil_type, "pF2"] + study_moisture <- meta$rel_moisture * meta$moisture_ref + object$f_time_norm <- f_time_norm_focus(meta$temperature, + moisture = study_moisture, field_moisture = field_moisture, + Q10 = Q10, walker = walker, f_na = f_na) + cat("$time_norm was set to\n") + print(object$f_time_norm) +} diff --git a/R/focus_soil_moisture.R b/R/focus_soil_moisture.R new file mode 100644 index 00000000..7b22fbcc --- /dev/null +++ b/R/focus_soil_moisture.R @@ -0,0 +1,11 @@ +utils::globalVariables("focus_soil_moisture") + +#' FOCUS default values for soil moisture contents at field capacity, MWHC and 1/3 bar +#' +#' @format A matrix with upper case USDA soil classes as row names, and water tension +#' ('pF1', 'pF2', 'pF 2.5') as column names +#' @source Anonymous (2014) Generic Guidance for Tier 1 FOCUS Ground Water Assessment +#' Version 2.2, May 2014 \url{https://esdac.jrc.ec.europa.eu/projects/ground-water} +#' @examples +#' focus_soil_moisture +"focus_soil_moisture" diff --git a/R/mkin_wide_to_long.R b/R/mkin_wide_to_long.R index 971f5273..3d0428d2 100644 --- a/R/mkin_wide_to_long.R +++ b/R/mkin_wide_to_long.R @@ -1,4 +1,4 @@ -if(getRversion() >= '2.15.1') utils::globalVariables(c("name", "time", "value")) +utils::globalVariables(c("name", "time", "value")) #' Convert a dataframe with observations over time into long format #' @@ -1,5 +1,5 @@ #' A dataset class for mkin -#' +#' #' @description #' At the moment this dataset class is hardly used in mkin. For example, #' mkinfit does not take mkinds datasets as argument, but works with dataframes @@ -7,12 +7,11 @@ #' provided by this package come as mkinds objects nevertheless. #' #' @importFrom R6 R6Class -#' @seealso The S3 printing method \code{\link{print.mkinds}} #' @examples -#' +#' #' mds <- mkinds$new("FOCUS A", FOCUS_2006_A) #' print(mds) -#' +#' #' @export mkinds <- R6Class("mkinds", public = list( @@ -62,15 +61,117 @@ mkinds <- R6Class("mkinds", ) #' Print mkinds objects -#' -#' @param x An \code{\link{mkinds}} object. +#' +#' @rdname mkinds +#' @param x An [mkinds] object. +#' @param data Should the data be printed? #' @param \dots Not used. #' @export -print.mkinds <- function(x, ...) { +print.mkinds <- function(x, data = FALSE, ...) { cat("<mkinds> with $title: ", x$title, "\n") cat("Observed compounds $observed: ", paste(x$observed, collapse = ", "), "\n") - cat("Sampling times $sampling_times: ", paste(x$sampling_times, collapse = ", "), "\n") + cat("Sampling times $sampling_times:\n") + cat(paste(x$sampling_times, collapse = ", "), "\n") cat("With a maximum of ", x$replicates, " replicates\n") if (!is.na(x$time_unit)) cat("Time unit: ", x$time_unit, "\n") if (!is.na(x$unit)) cat("Observation unit: ", x$unit, "\n") + if (data) print(mkin_long_to_wide(x$data)) +} + +#' A class for dataset groups for mkin +#' +#' @description +#' A container for working with datasets that share at least one compound, +#' so that combined evaluations are desirable. +#' +#' Time normalisation factors are initialised with a value of 1 for each +#' dataset if no data are supplied. +#' +#' @examples +#' +#' mdsg <- mkindsg$new("Experimental X", experimental_data_for_UBA_2019[6:10]) +#' print(mdsg) +#' print(mdsg, verbose = TRUE) +#' print(mdsg, verbose = TRUE, data = TRUE) +#' +#' @export +mkindsg <- R6Class("mkindsg", + public = list( + + #' @field title A title for the dataset group + title = NULL, + + #' @field ds A list of mkinds objects + ds = NULL, + + #' @field observed_n Occurrence counts of compounds in datasets + observed_n = NULL, + + #' @field f_time_norm Time normalisation factors + f_time_norm = NULL, + + #' @field meta A data frame with a row for each dataset, + #' containing additional information in the form + #' of categorical data (factors) or numerical data + #' (e.g. temperature, moisture, + #' or covariates like soil pH). + meta = NULL, + + #' @description + #' Create a new mkindsg object + #' @param title The title + #' @param ds A list of mkinds objects + #' @param f_time_norm Time normalisation factors + #' @param meta The meta data + initialize = function(title = "", ds, + f_time_norm = rep(1, length(ds)), meta) + { + self$title <- title + if (all(sapply(ds, inherits, "mkinds"))) { + self$ds <- ds + } else { + stop("Please supply a list of mkinds objects") + } + + all_observed <- unlist(lapply(ds, function(x) x$observed)) + observed <- factor(all_observed, levels = unique(all_observed)) + self$observed_n <- table(observed) + names(dimnames(self$observed_n)) <- NULL + self$f_time_norm <- f_time_norm + + if (!missing(meta)) { + self$meta <- meta + } + } + ) +) + +#' Print mkindsg objects +#' +#' @rdname mkindsg +#' @param x An [mkindsg] object. +#' @param verbose Should the mkinds objects be printed? +#' @param data Should the mkinds objects be printed with their data? +#' @param \dots Not used. +#' @export +print.mkindsg <- function(x, data = FALSE, verbose = data, ...) { + cat("<mkindsg> holding", length(x$ds), "mkinds objects\n") + cat("Title $title: ", x$title, "\n") + cat("Occurrene of observed compounds $observed_n:\n") + print(x$observed_n) + if (any(x$f_time_norm != 1)) { + cat("Time normalisation factors $f_time_norm:\n") + print(x$f_time_norm) + } + if (!is.null(x$meta)) { + cat("Meta information $meta:\n") + print(x$meta) + } + if (verbose) { + cat("\nDatasets $ds:") + for (ds in x$ds) { + cat("\n") + print(ds, data = data) + } + } } diff --git a/R/mkinerrmin.R b/R/mkinerrmin.R index f52692ba..388060d4 100644 --- a/R/mkinerrmin.R +++ b/R/mkinerrmin.R @@ -1,4 +1,4 @@ -if(getRversion() >= '2.15.1') utils::globalVariables(c("name", "value_mean")) +utils::globalVariables(c("name", "value_mean")) #' Calculate the minimum error to assume in order to pass the variance test #' diff --git a/R/mkinerrplot.R b/R/mkinerrplot.R index 36e22a43..02a12672 100644 --- a/R/mkinerrplot.R +++ b/R/mkinerrplot.R @@ -1,4 +1,4 @@ -if(getRversion() >= '2.15.1') utils::globalVariables(c("variable", "residual")) +utils::globalVariables(c("variable", "residual")) #' Function to plot squared residuals and the error model for an mkin object #' diff --git a/R/mkinfit.R b/R/mkinfit.R index 7fa1c56e..a6efc858 100644 --- a/R/mkinfit.R +++ b/R/mkinfit.R @@ -1,4 +1,4 @@ -if(getRversion() >= '2.15.1') utils::globalVariables(c("name", "time", "value")) +utils::globalVariables(c("name", "time", "value")) #' Fit a kinetic model to data with one or more state variables #' diff --git a/R/mkinresplot.R b/R/mkinresplot.R index 0bfdd02f..bad28ae8 100644 --- a/R/mkinresplot.R +++ b/R/mkinresplot.R @@ -1,4 +1,4 @@ -if(getRversion() >= '2.15.1') utils::globalVariables(c("variable", "residual")) +utils::globalVariables(c("variable", "residual")) #' Function to plot residuals stored in an mkin object #' diff --git a/R/plot.mixed.mmkin.R b/R/plot.mixed.mmkin.R index 903c3213..0d6b15c0 100644 --- a/R/plot.mixed.mmkin.R +++ b/R/plot.mixed.mmkin.R @@ -1,4 +1,4 @@ -if(getRversion() >= '2.15.1') utils::globalVariables("ds") +utils::globalVariables("ds") #' Plot predictions from a fitted nonlinear mixed model obtained via an mmkin row object #' diff --git a/R/plot.mkinfit.R b/R/plot.mkinfit.R index 48df9483..6ac4e45d 100644 --- a/R/plot.mkinfit.R +++ b/R/plot.mkinfit.R @@ -1,4 +1,4 @@ -if(getRversion() >= '2.15.1') utils::globalVariables(c("type", "variable", "observed")) +utils::globalVariables(c("type", "variable", "observed")) #' Plot the observed data and the fitted model of an mkinfit object #' @@ -102,7 +102,7 @@ saem.mmkin <- function(object, grDevices::png(tmp) } fit_time <- system.time({ - capture.output(f_saemix <- saemix::saemix(m_saemix, d_saemix, control), split = !quiet) + utils::capture.output(f_saemix <- saemix::saemix(m_saemix, d_saemix, control), split = !quiet) f_pred <- try(saemix::saemix.predict(f_saemix), silent = TRUE) if (!inherits(f_pred, "try-error")) { f_saemix <- f_pred |