qpcRo: qPCR

###############################################################################
#
# author: Christoph Kiefer
# email: christophak@bmb.sdu.dk
#
################################################################################
#'
#' Plot Cq values for quality control
#'
#' @param scheme A pipetting scheeme, e.g., generated by get.pipettingScheme.
#'     Generally a data.frame containing gene, cond, col and row information.
#' @param notPass Character vector containing wells, that are not to be included
#'     in the analysis.
#' @param rep_biol For which biological replicate should the Cq values
#'     be plotted.
#'
#' @import dplyr
#' @import ggplot2
#' @import purrr
#' @import ggpmisc
#' @import ggrepel
#'
#' @examples scheme <- get.pipettingScheme(c(paste0("gene_", 1:6), paste0("housekeeper_", 1:2)),
#'         c(0, 2, 4, "6_NTC", "6_CL", "6_iso", "7_NTC", "7_CL", "7_iso"), 2,
#'         c(2,2,2,3,3,3,3,3,3), 1)
#'     scheme <- import.LCcq("data-raw/example.txt", scheme)
#'     plot.cq(scheme, notPass = "N19")
#'
#' @export
plot.cq <- function(scheme, notPass = NULL, rep_biol = 1) {

    scheme %>%
        mutate(pos = paste0(row, col)) %>%
        group_by(gene, cond, repl_biol, repl_tech) %>%
        mutate(repl_pcr = as.factor(row_number())) %>%
        ungroup() %>%
        filter(!is.na(cq)) %>%
        mutate(pos = paste0(row, col)) %>%
        mutate(pass = if_else(pos %in% notPass, FALSE, TRUE)) %>%
        filter(repl_biol == rep_biol) %>%
        ggplot(aes(cond, -cq)) +
            geom_jitter(aes(colour = as.factor(repl_tech)), size = .5, width = .1, height = 0) +
            facet_wrap(~gene, scales = "free_y") +
        scale_color_discrete(guide = FALSE) +
        geom_text_repel(aes(label = pos, col = as.factor(if_else(pass, 2, 1))),
            max.iter = 50,
            family = "serif",
            segment.size = .25,
            size = 1.5)

}

#'
#' Take means from Cq values of pcr replicates
#'
#' @param scheme A pipetting scheeme, e.g., generated by get.pipettingScheme.
#'     Generally a data.frame containing gene, cond, col and row information.
#' @param notPass Character vector containing wells, that are not to be included
#'     in the analysis.
#'
#' @import dplyr
#'
#' @export
get.averageCq <- function(scheme, notPass = NULL) {
    scheme %>%
        mutate(pos = paste0(row, col)) %>%
        filter(!pos %in% notPass) %>%
        group_by(gene, cond, repl_biol, repl_tech) %>%
        summarise(cq = mean(cq, na.rm = TRUE)) %>%
        ungroup()

}

#'
#' Calculate delta Cq values
#'
#' @param scheme A pipetting scheeme, e.g., generated by get.pipettingScheme.
#'     Generally a data.frame containing gene, cond, col and row information.
#'
#' @import dplyr
#'
#' @export
get.dCq <- function(scheme, hkp) {
    ref <- scheme %>%
        filter(gene %in% hkp) %>%
        group_by(cond, repl_biol, repl_tech) %>%
        summarise(housekeeper = mean(cq)) %>%
        ungroup()

   scheme %>%
        filter(!gene %in% hkp) %>%
        left_join(ref, by = c("cond", "repl_biol", "repl_tech")) %>%
        mutate(dcq = cq - housekeeper)

}

#'
#' Take means from delta Cq values of technical, e.g. cell culture, replicates
#'
#' @param scheme A pipetting scheeme, e.g., generated by get.pipettingScheme.
#'     Generally a data.frame containing gene, cond, col and row information.
#'
#' @import dplyr
#' @import ggplot2
#'
#' @export
get.averageDCq <- function(scheme) {
    scheme %>%
        group_by(gene, cond, repl_biol) %>%
        summarise(dcq = mean(dcq, na.rm = TRUE))

}

#'
#' Take means from delta Cq values of technical, e.g. cell culture, replicates
#'
#' @param scheme A pipetting scheeme, e.g., generated by get.pipettingScheme.
#'     Generally a data.frame containing gene, cond, col and row information.
#'
#' @import dplyr
#'
#' @export
plot.dCq <- function(scheme) {
    scheme %>%
        ggplot(aes(cond, -dcq)) +
            geom_jitter(aes(colour = as.factor(repl_biol)), size = .5, width = .1, height = 0) +
            facet_wrap(~gene, scales = "free_y") +
        scale_color_discrete(guide = FALSE)

}


#'
#' Get delta Cq values
#'
#' @examples analyse.cq(scheme, "./data-raw/example.txt",
#'   hkp = c("housekeeper_1", "housekeeper_2"), silent = TRUE)
#'
#' @export
analyse.dcq <- function(scheme, file, hkp, output = ".", silent = FALSE,
    decimal_mark = '.') {

    df <- import.LCcq(file, scheme, decimal_mark = decimal_mark)

    if(silent == FALSE) {
        for (i in 1:max(df$repl_biol)) {
            plot.cq(df, rep_biol = i)
            ggsave(paste0(output, "/cq_rep_", i, ".pdf"), height = 10, width = 16.18)
        }
    }

    df <- df %>%
        get.averageCq() %>%
        get.dCq(hkp = hkp) %>%
        get.averageDCq()

    if(silent == FALSE) {
        write_csv(df, paste0(output, "/dcq.csv"))

        plot.dCq(df)
        ggsave(paste0(output, "/dcq.pdf"), height = 10, width = 16.18)
    }

    return(df)
}


#'
#' Calculate delta delta Cq values
#'
#' @export
get.ddcq <- function(scheme, reference = scheme$cond[1]) {
    scheme %>%
        group_by(gene, repl_biol) %>%
        mutate(ddcq = dcq - dcq[cond == reference]) %>%
        ungroup()
}

#'
#' Plot delta delta Cq values
#'
#' @export
plot.ddCq <- function(scheme, reorder = NULL) {

    if(!is.null(reorder)) {
        scheme <- scheme %>%
            mutate(cond = factor(cond, levels = reorder))
    }

    scheme %>%
        ggplot(aes(cond, -ddcq)) +
            geom_jitter(aes(colour = as.factor(repl_biol)), size = .5, width = .1, height = 0) +
            facet_wrap(~gene, scales = "free_y")

}

#'
#' Run the whole ddcq analysis
#'
#' @examples analyse.ddcq(scheme, "./data-raw/example.txt",
#'   hkp = c("housekeeper_1", "housekeeper_2"), silent = TRUE)
#'
#' @export
analyse.ddcq <- function(scheme, file, hkp, output = ".", silent = FALSE,
    reference = scheme$cond[1], decimal_mark = '.') {

    df <- scheme %>%
        analyse.dcq(file, hkp, output, silent, decimal_mark = decimal_mark) %>%
        get.ddcq()

    if(silent == FALSE) {
        plot.ddCq(df)
        ggsave(paste0(output, "/ddcq.pdf"), height = 10, width = 16.18)

        write_csv(df, paste0(output, "/ddcq.csv"))
    }

    return(df)

}

#'
#' Making nice breaks
#'
my_breaks <- function (n = 5, ...) {

    function(x) {
        dings <- .045 * (max(x) - min(x))
        extended_range_breaks_(min(x) + dings , max(x) - dings, n, ...)
    }
}

#'
#' Make your plot look nice
#'
#' @export
make.nice <- function(ggplot) {

    ggplot +
        theme_tufte() +
        geom_rangeframe(sides = 'l') +
        scale_y_continuous(breaks = my_breaks(),
            labels = function(x) round(x, 1),
            name = bquote('-'*Delta*Delta*'Cq')) +
        scale_color_viridis_d(guide = FALSE)
}