In onwhenrdy/pksimhelper: A PK-Sim/MoBi PBPK model evaluation package

knitr::knit_hooks$set(
   error = function(x, options) {
     paste('\n\n<div class="alert alert-danger">',
           gsub('##', '\n', gsub('^##\ Error', '**Error**', x)),
           '</div>', sep = '\n')
   },
   warning = function(x, options) {
     paste('\n\n<div class="alert alert-warning">',
           gsub('##', '\n', gsub('^##\ Warning:', '**Warning**', x)),
           '</div>', sep = '\n')
   },
   message = function(x, options) {
     paste('\n\n<div class="alert alert-info">',
           gsub('##', '\n', x),
           '</div>', sep = '\n')
   }
)

work.dir <- params$work.dir
PAGE.WITH <- as.numeric(params$max_pic_size)

# check for working dir
if (!file.exists(work.dir)) {
  knitr::knit_exit()
  stop(paste0("work.dir '", work.dir, "' does not exist"))
}

# check if at least one model was selected
model.dirs <- params$models
model.dirs <- model.dirs[!is.null(model.dirs) || !is.na(model.dirs)]
if (length(model.dirs) == 0) {
  knitr::knit_exit()
  stop("No model was set for evaluation")
}

if (length(model.dirs) > 4) {
  knitr::knit_exit()
  stop("The comparison of more than 4 models is not supported")
}

# check for model folders and pksh_eval.json
# all data that we need is now in model.data !!
model.dirs <- file.path(work.dir, model.dirs)
model.data <- list()
for (d in model.dirs) {
  if (!file.exists(d)) {
     knitr::knit_exit()
     stop(paste0("Model directory '", d, "' does not exist"))
  }

  meta.file <- file.path(d, "pksh_eval.json")
  if (!file.exists(meta.file)) {
    stop(paste0("pksh_eval.json in directory '", d, "' does not exist"))
  } else {
    meta.entry <- list(path = d,
                       meta = jsonlite::read_json(file.path(d, "pksh_eval.json")))
    model.data <- c(model.data, list(meta.entry))
  }
}

tabs <- if (params$use_tabs) "{.tabset}" else ""
PIC_SIZE <- floor(PAGE.WITH / length(model.data))

###########################################################################################
# Helper functions
suppressMessages(library(tidyverse))

plot_image <- function(file.name, model_name) {

  if (length(file.name) == 0 || !file.exists(file.name)) {
    plot.new()
    image(matrix(PIC_SIZE, PIC_SIZE), col = "white", axes = F)
    text(x = 0.5, y = 0, paste("Plot for", model_name, "is missing"), col = "red")
  }
  else { 
    cat(sprintf('![](%s){ width=%ipx }', URLencode(file.name), PIC_SIZE), sep = '\n\n')
    cat("  ")
  }  
}

get_model_names <- function(data) sapply(data, function(x) x$meta$model_name[[1]], simplify = T)
get_model_paths <- function(data) sapply(data, function(x) x$path, simplify = T)
print_comp_phrase <- function(names) cat(paste("Comparing __", paste(names, collapse = "__ vs __"), "__\n\n", sep = ""))

###########################################################################################

if (params$show_settings) {
  cat("# Settings")

  model_names <- get_model_names(model.data)
  model_paths <- get_model_paths(model.data)

  row.names <- c("Workding dir", paste("Model", 1:length(model.data))) 
  settings.table <- data.frame(keys = row.names,
                               names = c(" ", model_names),
                               paths = c(work.dir, model_paths))
  knitr::kable(settings.table, col.names = c("Key", "Name", "Path"))
}

if (params$show_individuals_lin || params$show_individuals_log) {
  cat("# Individual Time Profiles \n")

  model_names <- get_model_names(model.data)
  model_paths <- get_model_paths(model.data)
  profile_names <- model.data %>% purrr::map(function(x) names(x$meta$individuals)) %>% unlist() %>% unique()

  if (length(profile_names) == 0) {
    knitr::knit_exit()
    stop("No individual profiles found")
  }

  for (profile in profile_names) {
    cat(sprintf("## %s %s \n\n", profile, tabs))
    print_comp_phrase(model_names)

    if (params$show_individuals_lin) {
      cat("### Linear \n\n")
      for (i in 1:length(model.data)) {
        tmp_data <- model.data[[i]]$meta
        lin_file <- file.path(model_paths[i], tmp_data$individual[[profile]]$lin[[1]])
        plot_image(lin_file, model_names[i])
      }
      cat("<br/>\n\n")
    }

    if (params$show_individuals_log) {
      cat("### Logarithmic \n\n")
      for (i in 1:length(model.data)) {
        tmp_data <- model.data[[i]]$meta
        log_file <- file.path(model_paths[i], tmp_data$individual[[profile]]$log[[1]])
        plot_image(log_file, model_names[i])

      }
      cat("<br/>\n\n")
    }
  }
}

if (params$show_pred_obs) {
  cat(sprintf("# Prediction vs. Overved Plots %s \n\n", tabs))

  model_names <- get_model_names(model.data)
  model_paths <- get_model_paths(model.data)

  pred_obs_names <- model.data %>% purrr::map(function(x) names(x$meta$pred_obs)) %>% unlist() %>% unique()
  if (length(pred_obs_names) == 0) {
    knitr::knit_exit()
    stop("No prediction vs observed files found")
  }

   for (pred_obs in pred_obs_names) {
    cat(sprintf("## %s \n\n", pred_obs))
    print_comp_phrase(model_names)

    for (i in 1:length(model.data)) {
      tmp_data <- model.data[[i]]$meta
      file <- file.path(model_paths[i], tmp_data$pred_obs[[pred_obs]][[1]])
      plot_image(file, model_names[i])
    }
    cat("<br/>\n\n")
   }
}

if (params$show_gof) {
  cat(paste("# GoF Plots", tabs, "\n"))

  model_names <- get_model_names(model.data)
  model_paths <- get_model_paths(model.data)

  # cmax
  cat("## Cmax\n\n")
  print_comp_phrase(model_names)

  for (i in 1:length(model.data)) {
    tmp_data <- model.data[[i]]$meta
    file <- file.path(model_paths[i], tmp_data$gof$c_max[[1]])
    plot_image(file, model_names[i])
  }
  cat("<br/>\n\n")

  # AUC
  cat("## AUC\n\n")
  cat(paste("Comparing __", paste(model_names, collapse = "__ vs __"), "__\n\n", sep = ""))

  for (i in 1:length(model.data)) {
    tmp_data <- model.data[[i]]$meta
    file <- file.path(model_paths[i], tmp_data$gof$auc[[1]])
    plot_image(file, model_names[i])
  }
  cat("<br/>\n\n")
}

if (params$show_metrics) {
  cat(paste("# Metrics", tabs, "\n"))

  model_names <- get_model_names(model.data)
  model_paths <- get_model_paths(model.data)

  all_data <- list()
  for (i in 1:length(model.data)) {
    tmp_data <- model.data[[i]]$meta
    file <- file.path(model.data[[i]]$path, tmp_data$cor_stats[[1]])
    data <- read.table(file, header = T, sep = "\t")

    if (length(params$metrics) != 1 || tolower(params$metrics[1]) != "all") {
      data <- data[params$metrics]
    }

    tmp <- as.data.frame(t(data))
    all_data[[i]] <- tmp 
  }

  combined <- do.call("cbind", all_data)
  # remove a couple of not really interesting metrics
  cmax <- combined[seq(1, by = 3, length.out = length(model_names))]
  colnames(cmax) <- model_names

  auc <- combined[seq(3, by = 3, length.out = length(model_names))]
  colnames(auc) <- model_names
}

if (params$show_metrics) {
  cat("## Cmax\n\n")
  knitr::kable(cmax)
}

if (params$show_metrics) {
  cat("## AUC\n\n")
  knitr::kable(auc)
}

onwhenrdy/pksimhelper documentation built on Aug. 19, 2022, 7:56 a.m.