From 19c83ecd8ed5ff8f82ff37d33170abd6d9cb56c8 Mon Sep 17 00:00:00 2001 From: aduvermy <arnaud.duvermy@ens-lyon.fr> Date: Fri, 18 Aug 2023 17:38:55 +0200 Subject: [PATCH] update roc curve --- src/v4/HTRSIM/dev/flat_full_bis.Rmd | 126 +++++++++++++++++++++++++++- 1 file changed, 122 insertions(+), 4 deletions(-) diff --git a/src/v4/HTRSIM/dev/flat_full_bis.Rmd b/src/v4/HTRSIM/dev/flat_full_bis.Rmd index a6766ab..4c09449 100644 --- a/src/v4/HTRSIM/dev/flat_full_bis.Rmd +++ b/src/v4/HTRSIM/dev/flat_full_bis.Rmd @@ -5255,10 +5255,18 @@ roc_plot <- function(comparison_df, ...) { args <- lapply(list(...), function(x) if (!is.null(x)) ggplot2::sym(x)) #comparison_df$isDE <- factor(comparison_df$isDE, levels= c(TRUE, FALSE)) - ggplot2::ggplot(comparison_df, ggplot2::aes(d = !isDE , m = p.adj, !!!args )) + - plotROC::geom_roc(n.cuts = 0, labels = FALSE) + - ggplot2::theme_bw() + - ggplot2::ggtitle("ROC curve") + p <- ggplot2::ggplot(comparison_df, ggplot2::aes(d = !isDE , m = p.adj, !!!args )) + + plotROC::geom_roc(n.cuts = 0, labels = FALSE) + + ggplot2::theme_bw() + + ggplot2::ggtitle("ROC curve") + + ## -- annotation AUC + df_AUC <- plotROC::calc_auc(p)[c("from", "AUC")] + df_AUC$AUC <- round(df_AUC$AUC, digits = 3) + annotations <- do.call(paste, c(df_AUC, sep = " - AUC: ")) + annotations <- paste(annotations, collapse = "\n") + p <- p + ggplot2::annotate("text", x = .75, y = .25, label = annotations) + return(p) } @@ -5876,6 +5884,116 @@ test_that("wrapperDESeq2 function works correctly", { }) +``` + +## Anova + + +```{r function-anova, filename = "anova"} + +#' Handle ANOVA Errors +#' +#' This function handles ANOVA errors and warnings during the ANOVA calculation process. +#' +#' @param l_TMB A list of fitted glmmTMB models. +#' @param group A character string indicating the group for which ANOVA is calculated. +#' @param ... Additional arguments to be passed to the \code{car::Anova} function. +#' +#' @return A data frame containing ANOVA results for the specified group. +#' @export +#' +#' @examples +#' l_tmb <- fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length, +#' data = iris, group_by = "Species", n.cores = 1) +#' anova_res <- handleAnovaError(l_tmb, "setosa", type = "III") +#' +#' @importFrom car Anova +#' @export +handleAnovaError <- function(l_TMB, group, ...) { + tryCatch( + expr = { + withCallingHandlers( + car::Anova(l_TMB[[group]], ...), + warning = function(w) { + message(paste(Sys.time(), "warning for group", group, ":", conditionMessage(w))) + invokeRestart("muffleWarning") + }) + }, + error = function(e) { + message(paste(Sys.time(), "error for group", group, ":", conditionMessage(e))) + NULL + } + ) +} + + +#' Perform ANOVA on Multiple glmmTMB Models in Parallel +#' +#' This function performs analysis of variance (ANOVA) on a list of \code{glmmTMB} +#' models in parallel for different groups specified in the list. It returns a list +#' of ANOVA results for each group. +#' +#' @param l_glmmTMB A list of \code{glmmTMB} models, with model names corresponding to the groups. +#' @param ... Additional arguments passed to \code{\link[stats]{anova}} function. +#' +#' @return A list of ANOVA results for each group. +#' +#' @examples +#' # Perform ANOVA +#' data(iris) +#' l_tmb<- fitModelParallel( Sepal.Length ~ Sepal.Width + Petal.Length, +#' data = iris, group_by = "Species", n.cores = 1 ) +#' anov_res <- anovaParallel(l_tmb , type = "III") +#' @importFrom stats anova +#' @export +anovaParallel <- function(l_tmb, ...) { + l_group <- attributes(l_tmb)$names + lapply(setNames(l_group, l_group), function(group) handleAnovaError(l_tmb, group, ...)) +} + + +``` + + +```{r test-anova} + + +test_that("handleAnovaError return correct ouptut", { + data(iris) + l_tmb <- fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length, + data = iris, group_by = "Species", n.cores = 1) + anova_res <- handleAnovaError(l_tmb, "setosa", type = "III") + + expect_s3_class(anova_res, "data.frame") + expect_equal(nrow(anova_res), 3) # Number of levels +}) + +test_that("handleAnovaError return correct ouptut", { + data(iris) + l_tmb <- fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length, + data = iris, group_by = "Species", n.cores = 1) + anova_res <- handleAnovaError(l_tmb, "INALID_GROUP", type = "III") + + expect_null(anova_res) +}) + + + +test_that("anovaParallel returns valid ANOVA results", { + data(iris) + l_tmb <- fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length, + data = iris, group_by = "Species", n.cores = 1) + anov_res <- anovaParallel(l_tmb, type = "III") + + expect_is(anov_res, "list") + expect_equal(length(anov_res), length(unique(iris$Species))) + +}) + + + + + ``` -- GitLab