R/hypothesisTest.R

In prostarproteomics/DAPAR: Tools for the Differential Analysis of Proteins Abundance with R

#'
#' @title Density plots of logFC values
#' 
#' @description 
#' This function show the density plots of Fold Change (the same as calculated
#' by limma) for a list of the comparisons of conditions in a differential
#' analysis.
#'
#' @param df_logFC A dataframe that contains the logFC values
#'
#' @param threshold_LogFC The threshold on log(Fold Change) to
#' distinguish between differential and non-differential data
#'
#' @param pal xxx
#'
#' @return A highcharts density plot
#'
#' @author Samuel Wieczorek
#'
#' @examples
#' data(Exp1_R25_prot, package="DAPARdata")
#' obj <- Exp1_R25_prot[seq_len(100)]
#' level <- 'protein'
#' metacell.mask <- match.metacell(GetMetacell(obj), c("Missing POV", "Missing MEC"), level)
#' indices <- GetIndices_WholeMatrix(metacell.mask, op = ">=", th = 1)
#' obj <- MetaCellFiltering(obj, indices, cmd = "delete")
#' qData <- Biobase::exprs(obj$new)
#' sTab <- Biobase::pData(obj$new)
#' res <- limmaCompleteTest(qData, sTab, comp.type = "OnevsAll")
#' pal <- ExtendPalette(2, "Dark2")
#' hc_logFC_DensityPlot(res$logFC, threshold_LogFC = 1, pal = pal)
#'
#' @export
#' 
#' @import highcharter
#'
hc_logFC_DensityPlot <- function(df_logFC,
                                 threshold_LogFC = 0,
                                 pal = NULL) {
    
    pkgs.require(c("stats", "RColorBrewer", "grDevices"))
    
    if (threshold_LogFC < 0) {
        warning("The parameter 'threshold_LogFC' must be positive or equal 
            to zero.")
        return(NULL)
    }


    hc <- highcharter::highchart() %>%
        hc_title(text = "log(FC) repartition") %>%
        my_hc_chart(chartType = "spline", zoomType = "x") %>%
        hc_legend(enabled = TRUE) %>%
        hc_xAxis(
            title = list(text = "log(FC)"),
            plotBands = list(
                list(
                    from = -threshold_LogFC, 
                    to = threshold_LogFC, 
                    color = "lightgrey")
                ),
            plotLines = list(
                list(
                    color = "grey", 
                    width = 2, 
                    value = 0, 
                    zIndex = 5
                    )
                )
        ) %>%
        hc_yAxis(title = list(text = "Density")) %>%
        hc_tooltip(
            headerFormat = "",
            pointFormat = "<b> {series.name} </b>: {point.y} ",
            valueDecimals = 2
        ) %>%
        my_hc_ExportMenu(filename = "densityplot") %>%
        hc_plotOptions(
            series = list(
                animation = list(duration = 100),
                connectNulls = TRUE,
                marker = list(enabled = FALSE)
            )
        )

    if (is.null(df_logFC) || ncol(df_logFC) == 0) {
        return(hc)
    }

    myColors <- NULL
    if (is.null(pal)) {
        warning("Color palette set to default.")
        myColors <- ExtendPalette(ncol(df_logFC), "Paired")
    } else {
        if (length(pal) != ncol(df_logFC)) {
            warning("The color palette has not the same dimension as the 
                number of samples")
            myColors <- ExtendPalette(pal, "Paired")
        }
        myColors <- pal
    }

    nValues <- nrow(df_logFC) * ncol(df_logFC)
    nInf <- length(which(df_logFC <= -threshold_LogFC))
    nSup <- length(which(df_logFC >= threshold_LogFC))
    nInside <- length(which(abs(df_logFC) < threshold_LogFC))
    hc <- hc %>%
        hc_colors(myColors)

    maxY.inf <- NULL
    maxY.inside <- NULL
    maxY.sup <- NULL
    minX <- NULL
    maxX <- NULL


    for (i in seq_len(ncol(df_logFC))) {
        tmp <- stats::density(df_logFC[, i])
        ind <- tmp$y[which(tmp$x <= -threshold_LogFC)]
        maxY.inf <- max(maxY.inf, ifelse(length(ind) == 0, 0, ind))
        .ind1 <- which(tmp$x > -threshold_LogFC)
        .ind2 <- which(tmp$x < threshold_LogFC)
        maxY.inside <- max(maxY.inf, tmp$y[intersect(.ind1, .ind2)])
        ind <- tmp$y[which(tmp$x > threshold_LogFC)]
        maxY.sup <- max(
            maxY.sup, 
            ifelse(length(ind) == 0, tmp$y[length(tmp$y)], ind)
            )
        minX <- min(minX, tmp$x)
        maxX <- max(maxX, tmp$x)


        hc <- hc_add_series(hc,
            data.frame(x = tmp$x, y = tmp$y),
            name = colnames(df_logFC)[i]
        )
    }

    ## add annotations
    if (threshold_LogFC > 0) {
        hc <- hc %>% hc_add_annotation(
            labelOptions = list(
                shape = "connector",
                backgroundColor = "lightgrey",
                # verticalAlign = 'bottom',
                align = "left",
                # distance=0,
                style = list(
                    fontSize = "1.5em",
                    textOutline = "1px white"
                ),
                borderWidth = 0,
                x = 20
            ),
            labels = list(
                list(
                    point = list(
                        xAxis = 0,
                        yAxis = 0,
                        x = 0,
                        y = maxY.inside
                    ),
                    text = paste0("n Filtered out = ", 
                        nInside, "<br>(", 
                        round(100 * nInside / nValues, digits = 2), "%)")
                )
            )
        )
    }
    if (threshold_LogFC >= minX) {
        hc <- hc %>%
            hc_add_annotation(
                labelOptions = list(
                    shape = "connector",
                    backgroundColor = "rgba(255,255,255,0.5)",
                    verticalAlign = "top",
                    borderWidth = 0,
                    crop = TRUE,
                    style = list(
                        color = "blue",
                        fontSize = "1.5em",
                        textOutline = "1px white"
                    ),
                    y = -10
                ),
                labels = list(
                    list(
                        point = list(
                            xAxis = 0,
                            yAxis = 0,
                            x = mean(c(minX, -threshold_LogFC)),
                            y = maxY.inf
                        ),
                        text = paste0("nInf = ", nInf, "<br>(", 
                            round(100 * nInf / nValues, digits = 2), ")%")
                    )
                )
            )
    }

    if (threshold_LogFC <= maxX) {
        hc <- hc %>% hc_add_annotation(
            labelOptions = list(
                shape = "connector",
                backgroundColor = "blue",
                verticalAlign = "top",
                borderWidth = 0,
                style = list(
                    color = "blue",
                    fontSize = "1.5em",
                    textOutline = "1px white"
                ),
                y = -5
            ),
            labels = list(
                list(
                    point = list(
                        xAxis = 0,
                        yAxis = 0,
                        x = mean(c(maxX, threshold_LogFC)),
                        y = maxY.sup
                    ),
                    text = paste0("nSup = ", nSup, "<br>(", 
                        round(100 * nSup / nValues, digits = 2), ")%")
                )
            )
        )
    }



    return(hc)
}