diff options
Diffstat (limited to 'inst/rmarkdown/templates/hierarchical_kinetics')
3 files changed, 166 insertions, 0 deletions
diff --git a/inst/rmarkdown/templates/hierarchical_kinetics/skeleton/header.tex b/inst/rmarkdown/templates/hierarchical_kinetics/skeleton/header.tex new file mode 100644 index 00000000..a2b7ce83 --- /dev/null +++ b/inst/rmarkdown/templates/hierarchical_kinetics/skeleton/header.tex @@ -0,0 +1 @@ +\definecolor{shadecolor}{RGB}{248,248,248} diff --git a/inst/rmarkdown/templates/hierarchical_kinetics/skeleton/skeleton.Rmd b/inst/rmarkdown/templates/hierarchical_kinetics/skeleton/skeleton.Rmd new file mode 100644 index 00000000..74fae164 --- /dev/null +++ b/inst/rmarkdown/templates/hierarchical_kinetics/skeleton/skeleton.Rmd @@ -0,0 +1,162 @@ +--- +title: "Hierarchical kinetic modelling of degradation data" +author: +date: +output: mkin::hierarchical_kinetics +geometry: margin=2cm +--- + +```{r packages, cache = FALSE, message = FALSE} +library(mkin) +library(knitr) +library(saemix) +library(parallel) +library(readxl) +``` + +```{r n_cores, cache = FALSE, echo = FALSE} +n_cores <- detectCores() + +if (Sys.info()["sysname"] == "Windows") { + cl <- makePSOCKcluster(n_cores) +} else { + cl <- makeForkCluster(n_cores) +} +``` + +\clearpage + +# Introduction + +This report shows hierarchical kinetic modelling for ... +The data were obtained from ... + +```{r ds} +data_path <- system.file("testdata", "lambda-cyhalothrin_soil_efsa_2014.xlsx", package = "mkin") +ds <- read_spreadsheet(data_path, valid_datasets = c(1:4, 7:13)) +covariates <- attr(ds, "covariates") +``` + +The covariate data are shown below. + +```{r results = "asis", dependson = "ds", echo = FALSE} +kable(covariates, caption = "Covariate data for all datasets") +``` + +\clearpage + +The datasets with the residue time series are shown in the tables below. Please +refer to the spreadsheet for details like data sources, treatment of values +below reporting limits and time step normalisation factors. + +```{r results = "asis", dependson = "ds", echo = FALSE} +for (ds_name in names(ds)) { + print( + kable(mkin_long_to_wide(ds[[ds_name]]), + caption = paste("Dataset", ds_name), + booktabs = TRUE, row.names = FALSE)) + cat("\n\\clearpage\n") +} +``` + +# Parent only evaluations + +The following code performs separate fits of the candidate degradation models +to all datasets using constant variance and the two-component error model. + +```{r parent-sep, dependson = "ds"} +parent_deg_mods <- c("SFO", "FOMC", "DFOP", "SFORB") +errmods <- c(const = "constant variance", tc = "two-component error") +parent_sep_const <- mmkin( + parent_deg_mods, ds, + error_model = "const", + cluster = cl, quiet = TRUE) +parent_sep_tc <- update(parent_sep_const, error_model = "tc") +``` + +To select the parent model, the corresponding hierarchical fits are performed below. + +```{r parent-mhmkin, dependson = "parent-sep"} +parent_mhmkin <- mhmkin(list(parent_sep_const, parent_sep_tc), cluster = cl) +status(parent_mhmkin) |> kable() +``` + +All fits terminate without errors (status OK). The check for ill-defined +parameters shows that not all random effect parameters can be robustly +quantified. + +```{r dependson = "parent_mhmkin"} +illparms(parent_mhmkin) |> kable() +``` + +Therefore, the fits are updated, excluding random effects that were +ill-defined according to the `illparms` function. + +```{r parent-mhmkin-refined} +parent_mhmkin_refined <- update(parent_mhmkin, + no_random_effect = illparms(parent_mhmkin)) +status(parent_mhmkin_refined) |> kable() +``` + +The most suitable model is selected based on the AIC. + +```{r, dependson = "parent-mhmkin"} +aic_parent <- AIC(parent_mhmkin_refined) +min_aic <- which(aic_parent == min(aic_parent), arr.ind = TRUE) +best_degmod_parent <- rownames(aic_parent)[min_aic[1]] +best_errmod_parent <- colnames(aic_parent)[min_aic[2]] +anova(parent_mhmkin_refined) |> kable(digits = 1) +``` + +Based on the AIC, the combination of the `r best_degmod_parent` degradation +model with the error model `r errmods[best_errmod_parent]` is identified to +be most suitable for the degradation of the parent. The check below +confirms that no ill-defined parameters remain for this combined model. + +```{r dependson = "parent-mhmkin"} +illparms(parent_mhmkin_refined[[best_degmod_parent, best_errmod_parent]]) +``` + +The corresponding fit is shown below. + +```{r parent-best-full, dependson = "parent-mhmkin"} +plot(parent_mhmkin_refined[[best_degmod_parent, best_errmod_parent]]) +``` + +Detailed listings of the parent fits are shown in the Appendix. + +\clearpage + +# Appendix + +## Summaries of saem fits + +### Parent fits + +```{r listings-parent, results = "asis", echo = FALSE, dependson = "parent_mhmkin"} +for (deg_mod in parent_deg_mods) { + for (err_mod in c("const", "tc")) { + caption <- paste("Hierarchical", deg_mod, "fit with", errmods[err_mod]) + tex_listing(parent_mhmkin[[deg_mod, err_mod]], caption) + } +} +``` + +### Refined parent fits + +```{r listings-parent-refined, results = "asis", echo = FALSE, dependson = "parent_mhmkin_refined"} +for (deg_mod in parent_deg_mods) { + for (err_mod in c("const", "tc")) { + caption <- paste("Refined hierarchical", deg_mod, "fit with", errmods[err_mod]) + tex_listing(parent_mhmkin_refined[[deg_mod, err_mod]], caption) + } +} +``` + +## Session info + +```{r, echo = FALSE, cache = FALSE} +parallel::stopCluster(cl) +sessionInfo() +``` + diff --git a/inst/rmarkdown/templates/hierarchical_kinetics/template.yaml b/inst/rmarkdown/templates/hierarchical_kinetics/template.yaml new file mode 100644 index 00000000..d8ab6a4d --- /dev/null +++ b/inst/rmarkdown/templates/hierarchical_kinetics/template.yaml @@ -0,0 +1,3 @@ +name: Hierarchical kinetics +description: Hierarchical kinetic modelling of degradation data +create_dir: true |