exec/initial_drc_qc_database_setup.R

In hemoshear/assayr2: Tools for analyzing assay data

# Setup initial database for DRC QC warnings

library(magrittr)
library(dplyr)
library(ggplot2)
library(tibble)
library(RSQLite)
library(assayr2)


data_dir <- "/data/shiny/HST101_DRC"
print(file.path(data_dir, dir(data_dir)))
d <- purrr::map(file.path(data_dir, dir(data_dir, pattern = ".RDS")), 
                ~ readRDS(.))

#file.path(data_dir, dir(data_dir))

readRDS(file.path(data_dir, dir(data_dir))[1])

names(d) <- file.path(data_dir, dir(data_dir, pattern = ".RDS"))

filter_curves <- function(d) {
  curve_ids <- d[[1]]$curve_id %>%
    gsub("[0-9][0-9][0-9]\\_run[1-2]\\_", "", .) %>%
    gsub("\\-CoA", "", .) %>%
    unique
}

labeled <- purrr::map(d , ~ 
                        as.tibble(.) %>% 
                        tidyr::nest(., -targ, -tx_run) %>%
                        dplyr::filter(grepl("13", targ)))
names(labeled) <- names(d)
labeled <- do.call('rbind', labeled)
labeled %<>% dplyr::filter(grepl("Isobutyryl", targ))
labeled$data %<>% purrr::map(~ mutate(., tx_conc = as.numeric(tx_conc %>% as.character)))
labeled$data %<>% purrr::map( ~ mutate(., 
                                       tx_conc = ifelse(tx_conc == 0, 
                                                        new_zeros(tx_conc), tx_conc)))
prop_drc <- purrr::map(labeled$data, ~
                         drc::drm(conc_corrected ~ tx_conc, data = ., 
                                  fct = drc::LL.4(), 
                                  control = drc::drmc(errorm = F, useD=FALSE)))

names(prop_drc) <- c(labeled$tx_run)
warning_codes <- purrr::map(prop_drc, ~ ll4_qc_warnings(.))
warning_codes %>% unlist

w <- data.frame(tx_run = names(prop_drc),
                warning_code = warning_codes %>% unlist)

warning_db <- dbConnect(RSQLite::SQLite(), "/data/shiny/HST101_DRC/qc_database/qc-warnings.sqlite")
dbWriteTable(warning_db, "warning_codes", w, overwrite=T)

tx_run2warning_code <- function(tx_run) {
  if (length(tx_run) > 1) {
    tx_run_set <- purrr::map_chr(tx_run, ~ shQuote(.)) %>% paste(collapse=',')
  } else {
    tx_run_set <- shQuote(tx_run)
  }
  query <-  paste0("SELECT * FROM warning_codes WHERE tx_run IN (", tx_run_set, ")")
  dbGetQuery(warning_db, query)
}

# now create warning code lookup table
warning_message <- data.frame(
  warning_code = rep(1:9),
  warning_message =  c("unbounded CI for slope parameter",
                       "unbounded CI for lower asymptote",
                       "unbounded CI for upper asymptote",
                       "unbounded CI for IC50",
                       "outliers may be present", 
                       "outliers may be present",
                       "curve span << asymptote span",
                       "imprecise CI for IC50 estimate",
                       "noisy relative to historical")
)

dbWriteTable(warning_db, "warning_messages", warning_message, overwrite=T)
get_warning_message <- function(warning_codes) {
  warning_messages <- dbGetQuery(warning_db, "SELECT * FROM warning_messages")
  this_warning_code_set <- warning_codes %>% strsplit(split = "\\,") %>% .[[1]] %>% as.integer
  messages <- warning_messages[warning_messages$warning_code %in% this_warning_code_set,]$warning_message
  names(messages) <- warning_messages[warning_messages$warning_code %in% this_warning_code_set,]$warning_code
  messages
}

get_qc_warnings <- function(tx_run) {
  warning_codes <- tx_run2warning_code(tx_run)
  res <- purrr::map(warning_codes$warning_code, ~ get_warning_message(.))
  names(res) <- warning_codes$tx_run
  res
}

# get_qc_warnings(c('H004232_PAU0607', 'H003833_PAU0602-1'))