From 4d4c7adc28a23c8a122614775125e9fc7226bd4f Mon Sep 17 00:00:00 2001
From: aduvermy <arnaud.duvermy@ens-lyon.fr>
Date: Wed, 13 Mar 2024 10:45:08 +0100
Subject: [PATCH] relevelling ability
---
DESCRIPTION | 3 +-
NAMESPACE | 4 +
R/actual_interactionfixeffects.R | 2 +-
R/actual_mainfixeffects.R | 2 +-
R/anova.R | 2 +-
R/basal_expression_scaling.R | 4 +-
R/counts_plot.R | 2 +-
R/datafrommvrnorm_manipulations.R | 2 +-
R/datafromuser_manipulations.R | 2 +-
R/evaluate_dispersion.R | 2 +-
R/evaluation_identity.R | 6 +-
R/evaluation_withmixedeffect.R | 2 +-
R/fake-section-title.R | 2 +-
R/fitmodel.R | 13 +-
R/glance_glmmtmb.R | 2 +-
R/inferencetoexpected.R | 2 +-
R/mlmetrics.R | 2 +-
R/mock_rnaseq.R | 2 +-
R/plot_metrics.R | 2 +-
R/precision_recall.R | 2 +-
R/prepare_data2fit.R | 2 +-
R/receiver_operating_characteristic.R | 2 +-
R/rocr_functions.R | 2 +-
R/sequencing_depth_scaling.R | 2 +-
R/setcorrelation.R | 2 +-
R/simulation.R | 2 +-
R/simulation_initialization.R | 2 +-
R/simulation_report.R | 4 +-
R/subsetgenes.R | 2 +-
R/tidy_glmmtmb.R | 2 +-
R/update_fittedmodel.R | 180 ++++++++--
R/utils.R | 25 +-
R/waldtest.R | 2 +-
R/wrapper_dds.R | 4 +-
dev/flat_full.Rmd | 340 ++++++++++++++++--
man/detect_categoricals_vars.Rd | 32 ++
man/first_non_null_index.Rd | 21 ++
man/get_rsquare_2plot.Rd | 4 +-
man/relevel_list_tmb_frame.Rd | 37 ++
man/relevelling_factors.Rd | 32 ++
man/updateParallel.Rd | 25 +-
.../test-actual_interactionfixeffects.R | 2 +-
tests/testthat/test-actual_mainfixeffects.R | 2 +-
tests/testthat/test-anova.R | 2 +-
.../testthat/test-basal_expression_scaling.R | 2 +-
tests/testthat/test-counts_plot.R | 2 +-
.../test-datafrommvrnorm_manipulations.R | 2 +-
.../test-datafromuser_manipulations.R | 2 +-
tests/testthat/test-evaluate_dispersion.R | 2 +-
tests/testthat/test-evaluation_identity.R | 2 +-
.../test-evaluation_withmixedeffect.R | 2 +-
tests/testthat/test-fitmodel.R | 12 +-
tests/testthat/test-glance_glmmtmb.R | 2 +-
tests/testthat/test-mock_rnaseq.R | 2 +-
tests/testthat/test-plot_metrics.R | 2 +-
tests/testthat/test-precision_recall.R | 2 +-
tests/testthat/test-prepare_data2fit.R | 2 +-
.../test-receiver_operating_characteristic.R | 2 +-
.../testthat/test-sequencing_depth_scaling.R | 2 +-
tests/testthat/test-setcorrelation.R | 2 +-
tests/testthat/test-simulation.R | 2 +-
.../testthat/test-simulation_initialization.R | 2 +-
tests/testthat/test-simulation_report.R | 2 +-
tests/testthat/test-subsetgenes.R | 2 +-
tests/testthat/test-tidy_glmmtmb.R | 2 +-
tests/testthat/test-update_fittedmodel.R | 99 ++++-
tests/testthat/test-utils.R | 11 +-
tests/testthat/test-waldtest.R | 2 +-
tests/testthat/test-wrapper_dds.R | 3 +-
vignettes/03-rnaseq_analysis.Rmd | 12 +-
70 files changed, 826 insertions(+), 143 deletions(-)
create mode 100644 man/detect_categoricals_vars.Rd
create mode 100644 man/first_non_null_index.Rd
create mode 100644 man/relevel_list_tmb_frame.Rd
create mode 100644 man/relevelling_factors.Rd
diff --git a/DESCRIPTION b/DESCRIPTION
index c630ea9..a3b7659 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -40,6 +40,7 @@ Suggests:
testthat
VignetteBuilder:
knitr
+Config/fusen/version: 0.5.2
Encoding: UTF-8
Roxygen: list(markdown = TRUE)
-RoxygenNote: 7.2.2
+RoxygenNote: 7.3.1
diff --git a/NAMESPACE b/NAMESPACE
index 7430465..9761c66 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -32,6 +32,7 @@ export(convert2Factor)
export(correlation_matrix_2df)
export(countMatrix_2longDtf)
export(counts_plot)
+export(detect_categoricals_vars)
export(diagnostic_plot)
export(drop_randfx)
export(endsWithDigit)
@@ -46,6 +47,7 @@ export(fillInInteraction)
export(fillInVariable)
export(filter_dataframe)
export(findAttribute)
+export(first_non_null_index)
export(fitModel)
export(fitModelParallel)
export(fitUpdate)
@@ -135,6 +137,8 @@ export(prepare_dataParallel)
export(rbind_evaldata_tmb_dds)
export(rbind_model_params_and_dispersion)
export(recall)
+export(relevel_list_tmb_frame)
+export(relevelling_factors)
export(removeDigitsAtEnd)
export(removeDuplicatedWord)
export(renameColumns)
diff --git a/R/actual_interactionfixeffects.R b/R/actual_interactionfixeffects.R
index 81ce0a0..2498ba0 100644
--- a/R/actual_interactionfixeffects.R
+++ b/R/actual_interactionfixeffects.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Filter DataFrame
#'
diff --git a/R/actual_mainfixeffects.R b/R/actual_mainfixeffects.R
index 6642c14..eed250c 100644
--- a/R/actual_mainfixeffects.R
+++ b/R/actual_mainfixeffects.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Calculate average values by group
diff --git a/R/anova.R b/R/anova.R
index c2a2c75..63c4270 100644
--- a/R/anova.R
+++ b/R/anova.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Handle ANOVA Errors
diff --git a/R/basal_expression_scaling.R b/R/basal_expression_scaling.R
index c2ddc5d..83b7bfa 100644
--- a/R/basal_expression_scaling.R
+++ b/R/basal_expression_scaling.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
@@ -17,7 +17,7 @@
#' @examples
#' dtf <- data.frame(mu_ij = c(10, 20, 30, 15, 25, 35, 40, 5, 12, 22))
#' dtf_with_bins <- getBinExpression(dtf, n_bins = 3)
-#'
+#'
getBinExpression <- function(dtf_coef, n_bins){
col2bin <- "mu_ij"
bin_labels <- cut(dtf_coef[[col2bin]], n_bins, labels = paste("BinExpression", 1:n_bins, sep = "_"))
diff --git a/R/counts_plot.R b/R/counts_plot.R
index d61a6ef..1d19fa9 100644
--- a/R/counts_plot.R
+++ b/R/counts_plot.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Generate a density plot of gene counts
diff --git a/R/datafrommvrnorm_manipulations.R b/R/datafrommvrnorm_manipulations.R
index cbc046f..50f3b65 100644
--- a/R/datafrommvrnorm_manipulations.R
+++ b/R/datafrommvrnorm_manipulations.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' getInput2mvrnorm
#'
diff --git a/R/datafromuser_manipulations.R b/R/datafromuser_manipulations.R
index 1c3ee6d..d7c796f 100644
--- a/R/datafromuser_manipulations.R
+++ b/R/datafromuser_manipulations.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Get data from user
diff --git a/R/evaluate_dispersion.R b/R/evaluate_dispersion.R
index 0b30f99..cf7d4d5 100644
--- a/R/evaluate_dispersion.R
+++ b/R/evaluate_dispersion.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Get Dispersion Comparison
diff --git a/R/evaluation_identity.R b/R/evaluation_identity.R
index 76d9efa..fb74451 100644
--- a/R/evaluation_identity.R
+++ b/R/evaluation_identity.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
@@ -43,7 +43,9 @@ compute_rsquare <- function(data, grouping_by = c("from", "description") ){
#' @return A data frame with additional columns for labeling in the plot.
#' @export
#' @examples
-#' data_rsquare <- data.frame(from = c("A", "B", "C"), description = c("Desc1", "Desc2", "Desc3"), R2 = c(0.9, 0.8, 0.7))
+#' data_rsquare <- data.frame(from = c("A", "B", "C"),
+#' description = c("Desc1", "Desc2", "Desc3"),
+#' R2 = c(0.9, 0.8, 0.7))
#' result <- get_rsquare_2plot(data_rsquare)
get_rsquare_2plot <- function(data_rsquare){
data_rsquare$pos_x <- -Inf
diff --git a/R/evaluation_withmixedeffect.R b/R/evaluation_withmixedeffect.R
index 7a406bd..1264d43 100644
--- a/R/evaluation_withmixedeffect.R
+++ b/R/evaluation_withmixedeffect.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Check if the formula contains a mixed effect structure.
diff --git a/R/fake-section-title.R b/R/fake-section-title.R
index 72d65fe..4da9d9a 100644
--- a/R/fake-section-title.R
+++ b/R/fake-section-title.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' @name prediction-class
diff --git a/R/fitmodel.R b/R/fitmodel.R
index 4ae9cb7..01cd629 100644
--- a/R/fitmodel.R
+++ b/R/fitmodel.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Check if Data is Valid for Model Fitting
#'
@@ -115,8 +115,11 @@ is_fullrank <- function(metadata, formula) {
e <- eigen(crossprod(model_matrix), symmetric = TRUE, only.values = TRUE)$values
modelFullRank <- e[1] > 0 && abs(e[length(e)] / e[1]) > 1e-13
- if (!modelFullRank)
- stop("The model matrix is not full rank, so the model cannot be fit as specified. One or more variables or interaction terms in the design formula are linear combinations of the others and must be removed.")
+ if (!modelFullRank) {
+ warning("The model matrix is not full rank. One or more variables or interaction terms in the design formula are linear combinations of the others.")
+ return(FALSE)
+ }
+
return(TRUE)
}
@@ -135,6 +138,7 @@ is_fullrank <- function(metadata, formula) {
#' @examples
#' fitModel("mtcars" , formula = mpg ~ cyl + disp, data = mtcars)
fitModel <- function(group , formula, data, ...) {
+ is_fullrank(data, formula)
# Fit the model using glm.nb from the GLmmTMB package
model <- glmmTMB::glmmTMB(formula, ..., data = data )
@@ -242,7 +246,7 @@ parallel_fit <- function(groups, group_by, formula, data, n.cores = NULL,
l_data2parallel <- prepare_dataParallel(groups, group_by, data)
clust <- parallel::makeCluster(n.cores, outfile = log_file , type= cl_type )
- parallel::clusterExport(clust, c("fitModel"))
+ parallel::clusterExport(clust, c("fitModel", "is_fullrank", "drop_randfx"))
results_fit <- parallel::parLapply(clust, X = l_data2parallel,
fun = launchFit,
group_by = group_by, formula = formula, ...)
@@ -273,7 +277,6 @@ fitModelParallel <- function(formula, data, group_by, n.cores = NULL, cl_type =
## Some verification
isValidInput2fit(data, formula)
- is_fullrank(data, formula)
is_validGroupBy(data, group_by)
## -- print log location
diff --git a/R/glance_glmmtmb.R b/R/glance_glmmtmb.R
index 1fbea83..822fd42 100644
--- a/R/glance_glmmtmb.R
+++ b/R/glance_glmmtmb.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Extracts the summary statistics from a list of glmmTMB models.
diff --git a/R/inferencetoexpected.R b/R/inferencetoexpected.R
index b6128e5..e0ac974 100644
--- a/R/inferencetoexpected.R
+++ b/R/inferencetoexpected.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Compare the results of inference with the ground truth data.
diff --git a/R/mlmetrics.R b/R/mlmetrics.R
index dcf8a05..c2778d3 100644
--- a/R/mlmetrics.R
+++ b/R/mlmetrics.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/R/mock_rnaseq.R b/R/mock_rnaseq.R
index a2227cd..1ea3751 100644
--- a/R/mock_rnaseq.R
+++ b/R/mock_rnaseq.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Check the validity of the dispersion matrix
diff --git a/R/plot_metrics.R b/R/plot_metrics.R
index 7fd01f4..2b94c1e 100644
--- a/R/plot_metrics.R
+++ b/R/plot_metrics.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Subset the glance DataFrame based on selected variables.
diff --git a/R/precision_recall.R b/R/precision_recall.R
index c9ececc..b61bccd 100644
--- a/R/precision_recall.R
+++ b/R/precision_recall.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/R/prepare_data2fit.R b/R/prepare_data2fit.R
index 482c5b3..f173d2a 100644
--- a/R/prepare_data2fit.R
+++ b/R/prepare_data2fit.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Convert count matrix to long data frame
diff --git a/R/receiver_operating_characteristic.R b/R/receiver_operating_characteristic.R
index bb3754e..4dd2891 100644
--- a/R/receiver_operating_characteristic.R
+++ b/R/receiver_operating_characteristic.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/R/rocr_functions.R b/R/rocr_functions.R
index c9ececc..b61bccd 100644
--- a/R/rocr_functions.R
+++ b/R/rocr_functions.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/R/sequencing_depth_scaling.R b/R/sequencing_depth_scaling.R
index 7c0777b..f1bad83 100644
--- a/R/sequencing_depth_scaling.R
+++ b/R/sequencing_depth_scaling.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Scale Counts Table
diff --git a/R/setcorrelation.R b/R/setcorrelation.R
index 34a8377..c38666a 100644
--- a/R/setcorrelation.R
+++ b/R/setcorrelation.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Compute Covariation from Correlation and Standard Deviations
diff --git a/R/simulation.R b/R/simulation.R
index 4399ec4..6e25b6f 100644
--- a/R/simulation.R
+++ b/R/simulation.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Get input for simulation based on coefficients
#'
diff --git a/R/simulation_initialization.R b/R/simulation_initialization.R
index 65e22f6..416a0ae 100644
--- a/R/simulation_initialization.R
+++ b/R/simulation_initialization.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Initialize variable
#'
diff --git a/R/simulation_report.R b/R/simulation_report.R
index c3548f0..6955f61 100644
--- a/R/simulation_report.R
+++ b/R/simulation_report.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
@@ -85,7 +85,7 @@ is_truthLabels_valid <- function(eval_data, col_param = "description", col_truth
#' report <- evaluation_report(l_res, NULL, mock_data,
#' coeff_threshold = 0.67, alt_hypothesis = "greaterAbs")
#' }
-#'
+#'
evaluation_report <- function(list_tmb, dds, mock_obj, coeff_threshold, alt_hypothesis, alpha_risk = 0.05,
palette_color = c(DESeq2 = "#500472", HTRfit ="#79cbb8"),
palette_shape = c(DESeq2 = 17, HTRfit = 19),
diff --git a/R/subsetgenes.R b/R/subsetgenes.R
index 74171bc..351d045 100644
--- a/R/subsetgenes.R
+++ b/R/subsetgenes.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Subset Genes in Genomic Data
diff --git a/R/tidy_glmmtmb.R b/R/tidy_glmmtmb.R
index 7c2ebb8..435160d 100644
--- a/R/tidy_glmmtmb.R
+++ b/R/tidy_glmmtmb.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/R/update_fittedmodel.R b/R/update_fittedmodel.R
index 24ba2e8..4ce264c 100644
--- a/R/update_fittedmodel.R
+++ b/R/update_fittedmodel.R
@@ -1,21 +1,25 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Update glmmTMB models in parallel.
#'
-#' This function fits glmmTMB models in parallel using multiple cores, allowing for faster computation.
+#' This function updates glmmTMB models in parallel using multiple cores, allowing for faster computation.
+#' It updates the models with new reference labels if specified.
+#' It can also be used to fit a new formula or to change additional parameters of glmmTMB (param : "...").
#'
#' @param formula Formula for the GLMNB model.
#' @param list_tmb List of glmmTMB objects.
+#' @param reference_labels Vector of reference labels. Default is c(), selecting the first alphanumeric label as reference.
#' @param n.cores Number of cores to use for parallel processing. If NULL, the function will use all available cores.
#' @param log_file File path for the log output (default: Rtmpdir/htrfit.log).
-#' @param cl_type cluster type (defautl "PSOCK"). "FORK" is recommanded for linux.
+#' @param cl_type cluster type (default "PSOCK"). "FORK" is recommended for linux.
#' @param ... Additional arguments to be passed to the glmmTMB::glmmTMB function.
#' @export
#' @return A list of updated GLMNB models.
#'
#' @examples
+#' # -- Example usage: update formula
#' data(iris)
#' groups <- unique(iris$Species)
#' group_by <- "Species"
@@ -23,23 +27,154 @@
#' fitted_models <- fitModelParallel(formula, iris, group_by, n.cores = 1)
#' new_formula <- Sepal.Length ~ Sepal.Width
#' results <- updateParallel(new_formula, fitted_models, n.cores = 1)
-updateParallel <- function(formula, list_tmb, n.cores = NULL, cl_type = "PSOCK",
- log_file = paste(tempdir(check = FALSE), "htrfit.log", sep = "/"), ...) {
-
-
- isValidInput2fit(list_tmb[[1]]$frame, formula)
-
- is_fullrank(list_tmb[[1]]$frame, formula)
+#' # Example usage: update reference
+#' # -- Load the mtcars dataset
+#' data("mtcars")
+#' # -- Specify categorical variables
+#' mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+#' levels(mtcars$vs) <- c("V-shaped", "straight")
+#' mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+#' levels(mtcars$am) <- c("automatic", "manual")
+#' # -- For each group of number of cylinders:
+#' # -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+#' list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am,
+#' data = mtcars, group_by = "cyl", n.cores = 1)
+#' # -- Relevel transmission and engine shape variables
+#' list_tmb <- updateParallel(formula = mpg ~ hp + vs + am, list_tmb,
+#' reference_labels = c("straight", "manual"), n.cores = 1)
+updateParallel <- function (formula, list_tmb, reference_labels = c(), n.cores = NULL, cl_type = "PSOCK",
+ log_file = paste(tempdir(check = FALSE), "htrfit.log", sep = "/"),
+ ...) {
+ stopifnot(is.list(list_tmb))
+ non_null_idx <- first_non_null_index(list_tmb)
+ if (!is.null(non_null_idx)){
+ stopifnot(is.data.frame(list_tmb[[non_null_idx]]$frame))
+ isValidInput2fit(list_tmb[[non_null_idx]]$frame, formula)
+ list_tmb <- relevelling_factors(list_tmb, reference_labels)
+ message(paste("Log file location", log_file, sep = ": "))
+ list_tmb <- parallel_update(formula, list_tmb, n.cores, log_file, cl_type, ...)
+ clear_memory(except_obj = list_tmb)
+ }
+ return(list_tmb)
+}
+
+
+#' Re-levels categorical variables in a the frame of a list of glmmTMB objects.
+#'
+#' This function re-levels categorical variables in a list of glmmTMB objects using the specified reference labels.
+#'
+#' @param list_tmb List of glmmTMB objects.
+#' @param categorical_vars Names of the categorical variables to be re-leveled.
+#' @param ref_labels Vector of reference labels corresponding to the categorical variables.
+#' @return A list of glmmTMB objects with re-leveled categorical variables.
+#' @export
+#'
+#' @examples
+#' # Example usage:
+#' # -- Load the mtcars dataset
+#' data("mtcars")
+#' ## -- specify categorical var
+#' mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+#' levels(mtcars$vs) <- c("V-shaped", "straight")
+#' mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+#' levels(mtcars$am) <- c("automatic", "manual")
+#' # -- For each group of number of cylinders,
+#' # -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+#' list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am,
+#' data = mtcars, group_by = "cyl", n.cores = 1)
+#' # -- Relevel transmission and engine shape variables
+#' relevel_list_tmb_frame(list_tmb, c("am", "vs"), c("manual", "straight"))
+relevel_list_tmb_frame <- function(list_tmb, categorical_vars, ref_labels){
+ names(ref_labels) <- categorical_vars
+ lapply(list_tmb, function( tmb_obj ) {
+ if (!is.null(tmb_obj)){
+ for (categorical_var in names(ref_labels)) {
+ reference_label <- ref_labels[categorical_var]
+ tmb_obj$frame[[categorical_var]] <- stats::relevel(tmb_obj$frame[[categorical_var]], ref = reference_label)
+ }
+ }
+ return(tmb_obj)
+ })
+}
+
+
+#' Detects categorical variables based on reference labels in a glmmTMB object.
+#'
+#' This function detects categorical variables based on reference labels in a glmmTMB object's frame.
+#'
+#' @param tmb_frame The data frame of a glmmTMB object.
+#' @param ref_labels Vector of reference labels corresponding to categorical variables.
+#' @return Names of the categorical variables detected.
+#' @export
+#'
+#' @examples
+#' data("mtcars")
+#' ## -- specify categorical var
+#' mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+#' levels(mtcars$vs) <- c("V-shaped", "straight")
+#' mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+#' levels(mtcars$am) <- c("automatic", "manual")
+#' ## -- For each group of number of cylinder,
+#' ## -- explain fuel consumption with engine shape, Gross horsepower, and transmission type
+#' list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am,
+#' data = mtcars, group_by = "cyl" , n.cores = 1)
+#' detect_categoricals_vars(list_tmb[["6"]]$frame, c("straight", "manual"))
+detect_categoricals_vars <- function(tmb_frame, ref_labels){
+ idx_col <- c()
+ for (reference_label in ref_labels){
+
+ catego_var_idx_col <- unique(which(tmb_frame == reference_label, arr.ind = T)[, "col"])
+
+ if (length(catego_var_idx_col) > 1) {
+ message_err <- paste("Label", reference_label, " detected in the metadata across different columns.\nUnable to determine the correct columns for re-leveling.\nPlease ensure that reference labels are specific to individual columns for re-leveling")
+ stop(message_err)
+ }
+
+ if (length(catego_var_idx_col) == 0) {
+ ref_label_str <- paste(ref_labels , collapse = ", ")
+ message_err <- paste("Label", reference_label, "not found in metadata.")
+ stop(message_err)
+ }
+
+ idx_col <- c(idx_col, catego_var_idx_col)
+ }
- ## -- ## -- print log location
- message( paste("Log file location", log_file, sep =': ') )
-
- # Fit models update in parallel and capture the results
- results <- parallel_update(formula, list_tmb, n.cores, log_file, cl_type, ...)
- return(results)
+ return(colnames(tmb_frame)[idx_col])
+}
+
+#' Relevels factors in a list of glmmTMB objects using specified reference labels.
+#'
+#' This function re-levels factors in a list of glmmTMB objects based on the specified reference labels.
+#'
+#' @param list_tmb List of glmmTMB objects.
+#' @param reference_labels Vector of reference labels.
+#' @return A list of glmmTMB objects with re-leveled factors.
+#' @export
+#'
+#' @examples
+#' data("mtcars")
+#' ## -- specify categorical var
+#' mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+#' levels(mtcars$vs) <- c("V-shaped", "straight")
+#' mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+#' levels(mtcars$am) <- c("automatic", "manual")
+#' ## -- For each group of number of cylinder,
+#' ## -- explain fuel consumption with engine shape, Gross horsepower, and transmission type
+#' list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am,
+#' data = mtcars, group_by = "cyl", n.cores = 1 )
+#' relevelling_factors(list_tmb , c("straight", "manual"))
+relevelling_factors <- function(list_tmb, reference_labels ){
+ if (length(reference_labels) > 0){
+ l_categorical_vars <- detect_categoricals_vars(list_tmb[[1]]$frame, reference_labels)
+ list_tmb <- relevel_list_tmb_frame(list_tmb, l_categorical_vars, reference_labels)
+ }
+ return(list_tmb)
}
+
+
+
#' Internal function to fit glmmTMB models in parallel.
#'
#' This function is used internally by \code{\link{updateParallel}} to fit glmmTMB models in parallel.
@@ -63,18 +198,13 @@ updateParallel <- function(formula, list_tmb, n.cores = NULL, cl_type = "PSOCK",
parallel_update <- function(formula, list_tmb, n.cores = NULL,
log_file = paste(tempdir(check = FALSE), "htrfit.log", sep = "/"),
cl_type = "PSOCK" , ...) {
-
if (is.null(n.cores)) n.cores <- max(1, parallel::detectCores(logical = FALSE) - 1)
-
message(paste("CPU(s) number :", n.cores, sep = " "))
message(paste("Cluster type :", cl_type, sep = " "))
-
-
clust <- parallel::makeCluster(n.cores, type= cl_type, outfile = log_file)
- parallel::clusterExport(clust, c("launchUpdate", "fitUpdate"))
+ parallel::clusterExport(clust, c("launchUpdate", "fitUpdate", "is_fullrank", "drop_randfx"))
updated_res <- parallel::parLapply(clust, X = list_tmb, fun = launchUpdate , formula = formula, ...)
parallel::stopCluster(clust) ; invisible(gc(reset = T, verbose = F, full = T));
-
return(updated_res)
}
@@ -99,8 +229,8 @@ parallel_update <- function(formula, list_tmb, n.cores = NULL,
#' updated_model <- fitUpdate("setosa", fitted_models[[1]], new_formula)
fitUpdate <- function(group, glmm_obj, formula , ...){
data <- glmm_obj$frame
+ is_fullrank(data, formula)
resUpdt <- stats::update(glmm_obj, formula, ...)
-
resUpdt$frame <- data
## save groupID => avoid error in future update
resUpdt$groupId <- group
@@ -111,7 +241,6 @@ fitUpdate <- function(group, glmm_obj, formula , ...){
## control in ... => avoid error in future update
controlArgs <- additional_args[['control']]
if (!is.null(controlArgs)) resUpdt$call$control <- controlArgs
-
return(resUpdt)
}
@@ -135,11 +264,12 @@ fitUpdate <- function(group, glmm_obj, formula , ...){
#' new_formula <- Sepal.Length ~ Sepal.Width
#' updated_model <- launchUpdate(fitted_models[[1]], new_formula)
launchUpdate <- function(glmm_obj, formula, ...) {
+ if (is.null(glmm_obj)) return(NULL)
group <- glmm_obj$groupId
tryCatch(
expr = {
withCallingHandlers(
- fitUpdate(group ,glmm_obj, formula, ...),
+ fitUpdate(group, glmm_obj, formula, ...),
warning = function(w) {
message(paste(Sys.time(), "warning for group", group ,":", conditionMessage(w)))
invokeRestart("muffleWarning")
diff --git a/R/utils.R b/R/utils.R
index 8a9ef8b..e3e3bd3 100644
--- a/R/utils.R
+++ b/R/utils.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Join two data frames using data.table
#'
@@ -23,7 +23,25 @@ join_dtf <- function(d1, d2, k1, k2) {
return(dt_joined %>% as.data.frame())
}
-
+#' Finds the index of the first non-null element in a list.
+#'
+#' This function searches a list and returns the index of the first non-null element.
+#'
+#' @param lst The list to search.
+#' @return The index of the first non-null element, or NULL if no non-null element is found.
+#' @export
+#'
+#' @examples
+#' my_list <- list(NULL, NULL, 3, 5, NULL)
+#' first_non_null_index(my_list) # Returns 3
+first_non_null_index <- function(lst) {
+ for (i in seq_along(lst)) {
+ if (!is.null(lst[[i]])) {
+ return(i)
+ }
+ }
+ return(NULL)
+}
#' Clean Variable Name
#'
@@ -270,8 +288,7 @@ reorderColumns <- function(df, columnOrder) {
clear_memory <- function(except_obj){
- a <- rm(list = setdiff(ls(), except_obj)) ; gc( reset = TRUE, verbose = FALSE )
- return(invisible(NULL))
+ rm(list = setdiff(ls(), except_obj)) ; invisible(gc( reset = TRUE, verbose = FALSE ))
}
diff --git a/R/waldtest.R b/R/waldtest.R
index af7b2d0..e430a2f 100644
--- a/R/waldtest.R
+++ b/R/waldtest.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Wald test for hypothesis testing
diff --git a/R/wrapper_dds.R b/R/wrapper_dds.R
index 81b4f01..b2c050c 100644
--- a/R/wrapper_dds.R
+++ b/R/wrapper_dds.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
#' Wrapper Function for DESeq2 Analysis
@@ -32,7 +32,7 @@
#' result <- wrap_dds(dds, lfcThreshold = 1, altHypothesis = "greater")
#' @export
wrap_dds <- function(dds, lfcThreshold , altHypothesis, correction_method = "BH") {
- dds_full <- S4Vectors::mcols(dds) %>% as.data.frame()
+ dds_full <- as.data.frame(S4Vectors::mcols(dds))
## -- dispersion
message("INFO: The dispersion values from DESeq2 were reparametrized to their reciprocals (1/dispersion).")
diff --git a/dev/flat_full.Rmd b/dev/flat_full.Rmd
index 3ea6f7b..c99374a 100644
--- a/dev/flat_full.Rmd
+++ b/dev/flat_full.Rmd
@@ -48,7 +48,25 @@ join_dtf <- function(d1, d2, k1, k2) {
return(dt_joined %>% as.data.frame())
}
-
+#' Finds the index of the first non-null element in a list.
+#'
+#' This function searches a list and returns the index of the first non-null element.
+#'
+#' @param lst The list to search.
+#' @return The index of the first non-null element, or NULL if no non-null element is found.
+#' @export
+#'
+#' @examples
+#' my_list <- list(NULL, NULL, 3, 5, NULL)
+#' first_non_null_index(my_list) # Returns 3
+first_non_null_index <- function(lst) {
+ for (i in seq_along(lst)) {
+ if (!is.null(lst[[i]])) {
+ return(i)
+ }
+ }
+ return(NULL)
+}
#' Clean Variable Name
#'
@@ -295,8 +313,7 @@ reorderColumns <- function(df, columnOrder) {
clear_memory <- function(except_obj){
- a <- rm(list = setdiff(ls(), except_obj)) ; gc( reset = TRUE, verbose = FALSE )
- return(invisible(NULL))
+ rm(list = setdiff(ls(), except_obj)) ; invisible(gc( reset = TRUE, verbose = FALSE ))
}
@@ -304,6 +321,15 @@ clear_memory <- function(except_obj){
```{r tests-utils}
+
+# Test for first_non_null_index function
+test_that("first_non_null_index returns the correct index", {
+ lst <- list(NULL, NULL, 3, 5, NULL)
+ expect_equal(first_non_null_index(lst), 3)
+})
+
+
+
# Test unitaires pour la fonction join_dtf
test_that("join_dtf réalise la jointure correctement", {
# Création de données de test
@@ -2592,8 +2618,11 @@ is_fullrank <- function(metadata, formula) {
e <- eigen(crossprod(model_matrix), symmetric = TRUE, only.values = TRUE)$values
modelFullRank <- e[1] > 0 && abs(e[length(e)] / e[1]) > 1e-13
- if (!modelFullRank)
- stop("The model matrix is not full rank, so the model cannot be fit as specified. One or more variables or interaction terms in the design formula are linear combinations of the others and must be removed.")
+ if (!modelFullRank) {
+ warning("The model matrix is not full rank. One or more variables or interaction terms in the design formula are linear combinations of the others.")
+ return(FALSE)
+ }
+
return(TRUE)
}
@@ -2612,6 +2641,7 @@ is_fullrank <- function(metadata, formula) {
#' @examples
#' fitModel("mtcars" , formula = mpg ~ cyl + disp, data = mtcars)
fitModel <- function(group , formula, data, ...) {
+ is_fullrank(data, formula)
# Fit the model using glm.nb from the GLmmTMB package
model <- glmmTMB::glmmTMB(formula, ..., data = data )
@@ -2719,7 +2749,7 @@ parallel_fit <- function(groups, group_by, formula, data, n.cores = NULL,
l_data2parallel <- prepare_dataParallel(groups, group_by, data)
clust <- parallel::makeCluster(n.cores, outfile = log_file , type= cl_type )
- parallel::clusterExport(clust, c("fitModel"))
+ parallel::clusterExport(clust, c("fitModel", "is_fullrank", "drop_randfx"))
results_fit <- parallel::parLapply(clust, X = l_data2parallel,
fun = launchFit,
group_by = group_by, formula = formula, ...)
@@ -2750,7 +2780,6 @@ fitModelParallel <- function(formula, data, group_by, n.cores = NULL, cl_type =
## Some verification
isValidInput2fit(data, formula)
- is_fullrank(data, formula)
is_validGroupBy(data, group_by)
## -- print log location
@@ -2858,8 +2887,9 @@ test_that("Detect rank-deficient model matrix and throw error", {
z = factor(rep(c("zA","zB"), each = 5)),
y = rnorm(10))
formula <- y ~ x + w + z + y:w
- expect_error(is_fullrank(metadata, formula),
- regexp = "The model matrix is not full rank, so the model cannot be fit as specified.")
+ expect_warning(is_fullrank(metadata, formula))
+ res <- suppressWarnings(is_fullrank(metadata, formula))
+ expect_false(res)
})
# Test if a rank-deficient model matrix is detected and throws an error
@@ -2869,8 +2899,9 @@ test_that("Detect rank-deficient model matrix and throw error (with random eff)"
z = factor(rep(c("zA","zB"), each = 5)),
y = rnorm(10))
formula <- y ~ x + w + z + y:w + (1 | w)
- expect_error(is_fullrank(metadata, formula),
- regexp = "The model matrix is not full rank, so the model cannot be fit as specified.")
+ expect_warning(is_fullrank(metadata, formula))
+ res <- suppressWarnings(is_fullrank(metadata, formula))
+ expect_false(res)
})
# Test if a rank-deficient model matrix is detected and throws an error
@@ -3003,18 +3034,22 @@ test_that("fitModelParallel fits models in parallel for each group and returns a
#' Update glmmTMB models in parallel.
#'
-#' This function fits glmmTMB models in parallel using multiple cores, allowing for faster computation.
+#' This function updates glmmTMB models in parallel using multiple cores, allowing for faster computation.
+#' It updates the models with new reference labels if specified.
+#' It can also be used to fit a new formula or to change additional parameters of glmmTMB (param : "...").
#'
#' @param formula Formula for the GLMNB model.
#' @param list_tmb List of glmmTMB objects.
+#' @param reference_labels Vector of reference labels. Default is c(), selecting the first alphanumeric label as reference.
#' @param n.cores Number of cores to use for parallel processing. If NULL, the function will use all available cores.
#' @param log_file File path for the log output (default: Rtmpdir/htrfit.log).
-#' @param cl_type cluster type (defautl "PSOCK"). "FORK" is recommanded for linux.
+#' @param cl_type cluster type (default "PSOCK"). "FORK" is recommended for linux.
#' @param ... Additional arguments to be passed to the glmmTMB::glmmTMB function.
#' @export
#' @return A list of updated GLMNB models.
#'
#' @examples
+#' # -- Example usage: update formula
#' data(iris)
#' groups <- unique(iris$Species)
#' group_by <- "Species"
@@ -3022,22 +3057,153 @@ test_that("fitModelParallel fits models in parallel for each group and returns a
#' fitted_models <- fitModelParallel(formula, iris, group_by, n.cores = 1)
#' new_formula <- Sepal.Length ~ Sepal.Width
#' results <- updateParallel(new_formula, fitted_models, n.cores = 1)
-updateParallel <- function(formula, list_tmb, n.cores = NULL, cl_type = "PSOCK",
- log_file = paste(tempdir(check = FALSE), "htrfit.log", sep = "/"), ...) {
-
-
- isValidInput2fit(list_tmb[[1]]$frame, formula)
-
- is_fullrank(list_tmb[[1]]$frame, formula)
+#' #' # Example usage: update reference
+#' # -- Load the mtcars dataset
+#' data("mtcars")
+#' # -- Specify categorical variables
+#' mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+#' levels(mtcars$vs) <- c("V-shaped", "straight")
+#' mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+#' levels(mtcars$am) <- c("automatic", "manual")
+#' # -- For each group of number of cylinders:
+#' # -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+#' list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am,
+#' data = mtcars, group_by = "cyl", n.cores = 1)
+#' # -- Relevel transmission and engine shape variables
+#' list_tmb <- updateParallel(formula = mpg ~ hp + vs + am, list_tmb,
+#' reference_labels = c("straight", "manual"), n.cores = 1)
+updateParallel <- function (formula, list_tmb, reference_labels = c(), n.cores = NULL, cl_type = "PSOCK",
+ log_file = paste(tempdir(check = FALSE), "htrfit.log", sep = "/"),
+ ...) {
+ stopifnot(is.list(list_tmb))
+ non_null_idx <- first_non_null_index(list_tmb)
+ if (!is.null(non_null_idx)){
+ stopifnot(is.data.frame(list_tmb[[non_null_idx]]$frame))
+ isValidInput2fit(list_tmb[[non_null_idx]]$frame, formula)
+ list_tmb <- relevelling_factors(list_tmb, reference_labels)
+ message(paste("Log file location", log_file, sep = ": "))
+ list_tmb <- parallel_update(formula, list_tmb, n.cores, log_file, cl_type, ...)
+ clear_memory(except_obj = list_tmb)
+ }
+ return(list_tmb)
+}
+
+
+#' Re-levels categorical variables in a the frame of a list of glmmTMB objects.
+#'
+#' This function re-levels categorical variables in a list of glmmTMB objects using the specified reference labels.
+#'
+#' @param list_tmb List of glmmTMB objects.
+#' @param categorical_vars Names of the categorical variables to be re-leveled.
+#' @param ref_labels Vector of reference labels corresponding to the categorical variables.
+#' @return A list of glmmTMB objects with re-leveled categorical variables.
+#' @export
+#'
+#' @examples
+#' # Example usage:
+#' # -- Load the mtcars dataset
+#' data("mtcars")
+#' ## -- specify categorical var
+#' mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+#' levels(mtcars$vs) <- c("V-shaped", "straight")
+#' mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+#' levels(mtcars$am) <- c("automatic", "manual")
+#' # -- For each group of number of cylinders,
+#' # -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+#' list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am,
+#' data = mtcars, group_by = "cyl", n.cores = 1)
+#' # -- Relevel transmission and engine shape variables
+#' relevel_list_tmb_frame(list_tmb, c("am", "vs"), c("manual", "straight"))
+relevel_list_tmb_frame <- function(list_tmb, categorical_vars, ref_labels){
+ names(ref_labels) <- categorical_vars
+ lapply(list_tmb, function( tmb_obj ) {
+ if (!is.null(tmb_obj)){
+ for (categorical_var in names(ref_labels)) {
+ reference_label <- ref_labels[categorical_var]
+ tmb_obj$frame[[categorical_var]] <- stats::relevel(tmb_obj$frame[[categorical_var]], ref = reference_label)
+ }
+ }
+ return(tmb_obj)
+ })
+}
+
+
+#' Detects categorical variables based on reference labels in a glmmTMB object.
+#'
+#' This function detects categorical variables based on reference labels in a glmmTMB object's frame.
+#'
+#' @param tmb_frame The data frame of a glmmTMB object.
+#' @param ref_labels Vector of reference labels corresponding to categorical variables.
+#' @return Names of the categorical variables detected.
+#' @export
+#'
+#' @examples
+#' data("mtcars")
+#' ## -- specify categorical var
+#' mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+#' levels(mtcars$vs) <- c("V-shaped", "straight")
+#' mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+#' levels(mtcars$am) <- c("automatic", "manual")
+#' ## -- For each group of number of cylinder,
+#' ## -- explain fuel consumption with engine shape, Gross horsepower, and transmission type
+#' list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am,
+#' data = mtcars, group_by = "cyl" , n.cores = 1)
+#' detect_categoricals_vars(list_tmb[["6"]]$frame, c("straight", "manual"))
+detect_categoricals_vars <- function(tmb_frame, ref_labels){
+ idx_col <- c()
+ for (reference_label in ref_labels){
+
+ catego_var_idx_col <- unique(which(tmb_frame == reference_label, arr.ind = T)[, "col"])
+
+ if (length(catego_var_idx_col) > 1) {
+ message_err <- paste("Label", reference_label, " detected in the metadata across different columns.\nUnable to determine the correct columns for re-leveling.\nPlease ensure that reference labels are specific to individual columns for re-leveling")
+ stop(message_err)
+ }
+
+ if (length(catego_var_idx_col) == 0) {
+ ref_label_str <- paste(ref_labels , collapse = ", ")
+ message_err <- paste("Label", reference_label, "not found in metadata.")
+ stop(message_err)
+ }
+
+ idx_col <- c(idx_col, catego_var_idx_col)
+ }
- ## -- ## -- print log location
- message( paste("Log file location", log_file, sep =': ') )
-
- # Fit models update in parallel and capture the results
- results <- parallel_update(formula, list_tmb, n.cores, log_file, cl_type, ...)
- return(results)
+ return(colnames(tmb_frame)[idx_col])
}
+#' Relevels factors in a list of glmmTMB objects using specified reference labels.
+#'
+#' This function re-levels factors in a list of glmmTMB objects based on the specified reference labels.
+#'
+#' @param list_tmb List of glmmTMB objects.
+#' @param reference_labels Vector of reference labels.
+#' @return A list of glmmTMB objects with re-leveled factors.
+#' @export
+#'
+#' @examples
+#' data("mtcars")
+#' ## -- specify categorical var
+#' mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+#' levels(mtcars$vs) <- c("V-shaped", "straight")
+#' mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+#' levels(mtcars$am) <- c("automatic", "manual")
+#' ## -- For each group of number of cylinder,
+#' ## -- explain fuel consumption with engine shape, Gross horsepower, and transmission type
+#' list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am,
+#' data = mtcars, group_by = "cyl", n.cores = 1 )
+#' relevelling_factors(list_tmb , c("straight", "manual"))
+relevelling_factors <- function(list_tmb, reference_labels ){
+ if (length(reference_labels) > 0){
+ l_categorical_vars <- detect_categoricals_vars(list_tmb[[1]]$frame, reference_labels)
+ list_tmb <- relevel_list_tmb_frame(list_tmb, l_categorical_vars, reference_labels)
+ }
+ return(list_tmb)
+}
+
+
+
+
#' Internal function to fit glmmTMB models in parallel.
#'
@@ -3062,18 +3228,13 @@ updateParallel <- function(formula, list_tmb, n.cores = NULL, cl_type = "PSOCK",
parallel_update <- function(formula, list_tmb, n.cores = NULL,
log_file = paste(tempdir(check = FALSE), "htrfit.log", sep = "/"),
cl_type = "PSOCK" , ...) {
-
if (is.null(n.cores)) n.cores <- max(1, parallel::detectCores(logical = FALSE) - 1)
-
message(paste("CPU(s) number :", n.cores, sep = " "))
message(paste("Cluster type :", cl_type, sep = " "))
-
-
clust <- parallel::makeCluster(n.cores, type= cl_type, outfile = log_file)
- parallel::clusterExport(clust, c("launchUpdate", "fitUpdate"))
+ parallel::clusterExport(clust, c("launchUpdate", "fitUpdate", "is_fullrank", "drop_randfx"))
updated_res <- parallel::parLapply(clust, X = list_tmb, fun = launchUpdate , formula = formula, ...)
parallel::stopCluster(clust) ; invisible(gc(reset = T, verbose = F, full = T));
-
return(updated_res)
}
@@ -3098,8 +3259,8 @@ parallel_update <- function(formula, list_tmb, n.cores = NULL,
#' updated_model <- fitUpdate("setosa", fitted_models[[1]], new_formula)
fitUpdate <- function(group, glmm_obj, formula , ...){
data <- glmm_obj$frame
+ is_fullrank(data, formula)
resUpdt <- stats::update(glmm_obj, formula, ...)
-
resUpdt$frame <- data
## save groupID => avoid error in future update
resUpdt$groupId <- group
@@ -3110,7 +3271,6 @@ fitUpdate <- function(group, glmm_obj, formula , ...){
## control in ... => avoid error in future update
controlArgs <- additional_args[['control']]
if (!is.null(controlArgs)) resUpdt$call$control <- controlArgs
-
return(resUpdt)
}
@@ -3134,11 +3294,12 @@ fitUpdate <- function(group, glmm_obj, formula , ...){
#' new_formula <- Sepal.Length ~ Sepal.Width
#' updated_model <- launchUpdate(fitted_models[[1]], new_formula)
launchUpdate <- function(glmm_obj, formula, ...) {
+ if (is.null(glmm_obj)) return(NULL)
group <- glmm_obj$groupId
tryCatch(
expr = {
withCallingHandlers(
- fitUpdate(group ,glmm_obj, formula, ...),
+ fitUpdate(group, glmm_obj, formula, ...),
warning = function(w) {
message(paste(Sys.time(), "warning for group", group ,":", conditionMessage(w)))
invokeRestart("muffleWarning")
@@ -3197,8 +3358,105 @@ test_that("updateParallel function returns correct results", {
n.cores = 1,
family = glmmTMB::ziGamma(link = "inverse")))
expect_s3_class(updated_updated_model$setosa$call$family, "family")
+
+
+ ## -- update label reference
+ data("mtcars")
+ # -- Specify categorical variables
+ mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+ levels(mtcars$vs) <- c("V-shaped", "straight")
+ mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+ levels(mtcars$am) <- c("automatic", "manual")
+ # -- For each group of number of cylinders:
+ # -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+ list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am, data = mtcars, group_by = "cyl", n.cores = 1)
+ # -- Relevel transmission and engine shape variables
+ reference_labels <- c("straight", "manual")
+ result <- updateParallel(formula = mpg ~ hp + vs + am, list_tmb = list_tmb, reference_labels = reference_labels, n.cores = 1)
+
+ # Check if the returned list has the same length as the mock list
+ expect_equal(length(result), length(list_tmb))
+ expect_equal(levels(result[["6"]]$frame$am)[1] , "manual")
+ expect_equal(levels(result[["4"]]$frame$am)[1] , "manual")
+ expect_equal(levels(result[["6"]]$frame$vs)[1] , "straight")
+ expect_equal(levels(result[["4"]]$frame$vs)[1] , "straight")
+
+})
+
+
+
+# Test for relevel_list_tmb_frame function
+test_that("relevel_list_tmb_frame re-levels categorical variables correctly", {
+ data("mtcars")
+ # -- Specify categorical variables
+ mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+ levels(mtcars$vs) <- c("V-shaped", "straight")
+ mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+ levels(mtcars$am) <- c("automatic", "manual")
+ # -- For each group of number of cylinders:
+ # -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+ list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am, data = mtcars, group_by = "cyl", n.cores = 1)
+ # -- Relevel transmission and engine shape variables
+ reference_labels <- c("straight", "manual")
+ result <- relevel_list_tmb_frame(list_tmb, c("vs", "am"), reference_labels)
+
+ # Check if the returned list has the same length as the mock list
+ expect_equal(length(result), length(list_tmb))
+ # Check if categorical variables have been re-leveled correctly
+ expect_equal(levels(result[["6"]]$frame$am)[1] , "manual")
+ expect_equal(levels(result[["4"]]$frame$am)[1] , "manual")
+ expect_equal(levels(result[["6"]]$frame$vs)[1] , "straight")
+ expect_equal(levels(result[["4"]]$frame$vs)[1] , "straight")
+})
+
+# Test for relevelling_factors function
+test_that("relevelling_factors relevel categorical variables correctly", {
+ data("mtcars")
+ # -- Specify categorical variables
+ mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+ levels(mtcars$vs) <- c("V-shaped", "straight")
+ mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+ levels(mtcars$am) <- c("automatic", "manual")
+ # -- For each group of number of cylinders:
+ # -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+ list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am, data = mtcars, group_by = "cyl", n.cores = 1)
+ # -- Relevel transmission and engine shape variables
+ reference_labels <- c("straight", "manual")
+ result <- relevelling_factors(list_tmb, reference_labels)
+
+ # Check if the returned list has the same length as the mock list
+ expect_equal(length(result), length(list_tmb))
+ # Check if categorical variables have been re-leveled correctly
+ expect_equal(levels(result[["6"]]$frame$am)[1] , "manual")
+ expect_equal(levels(result[["4"]]$frame$am)[1] , "manual")
+ expect_equal(levels(result[["6"]]$frame$vs)[1] , "straight")
+ expect_equal(levels(result[["4"]]$frame$vs)[1] , "straight")
+})
+
+
+
+# Test for detect_categoricals_vars function
+test_that("detect_categoricals_vars detect categorical variables correctly", {
+ data("mtcars")
+ # -- Specify categorical variables
+ mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+ levels(mtcars$vs) <- c("V-shaped", "straight")
+ mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+ levels(mtcars$am) <- c("automatic", "manual")
+ # -- For each group of number of cylinders:
+ # -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+ list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am, data = mtcars, group_by = "cyl", n.cores = 1)
+ # -- Relevel transmission and engine shape variables
+ reference_labels <- c("straight", "manual")
+ result <- detect_categoricals_vars(list_tmb[["6"]]$frame, reference_labels)
+
+ # Check if the returned list has the same length as the mock list
+ expect_equal(length(result), length(reference_labels))
+ # Check if categorical variables have been re-leveled correctly
+ expect_equal(result, c("vs", "am"))
})
+
# Test parallel_update function
test_that("parallel_update function returns correct results", {
# Load the required data
@@ -6011,7 +6269,9 @@ compute_rsquare <- function(data, grouping_by = c("from", "description") ){
#' @return A data frame with additional columns for labeling in the plot.
#' @export
#' @examples
-#' data_rsquare <- data.frame(from = c("A", "B", "C"), description = c("Desc1", "Desc2", "Desc3"), R2 = c(0.9, 0.8, 0.7))
+#' data_rsquare <- data.frame(from = c("A", "B", "C"),
+#' description = c("Desc1", "Desc2", "Desc3"),
+#' R2 = c(0.9, 0.8, 0.7))
#' result <- get_rsquare_2plot(data_rsquare)
get_rsquare_2plot <- function(data_rsquare){
data_rsquare$pos_x <- -Inf
@@ -8405,7 +8665,7 @@ test_that("get_eval_data returns correct output", {
#' result <- wrap_dds(dds, lfcThreshold = 1, altHypothesis = "greater")
#' @export
wrap_dds <- function(dds, lfcThreshold , altHypothesis, correction_method = "BH") {
- dds_full <- S4Vectors::mcols(dds) %>% as.data.frame()
+ dds_full <- as.data.frame(S4Vectors::mcols(dds))
## -- dispersion
message("INFO: The dispersion values from DESeq2 were reparametrized to their reciprocals (1/dispersion).")
@@ -8653,6 +8913,7 @@ test_that("wrapperDESeq2 function works correctly", {
})
+
```
@@ -9342,7 +9603,7 @@ test_that("calculate_actualMixed calculates actual mixed effects as expected", {
```
```{r development-inflate, eval=FALSE}
-setwd("/Users/ex_dya/Documents/LBMC/HTRfit/")
+setwd("/home/adminarnaud/Documents/HTRfit/")
#usethis::create_package(path = "/Users/ex_dya/Documents/LBMC/HTRfit/")
fusen::fill_description(fields = list(Title = "HTRfit"), overwrite = T)
usethis::use_gpl_license(version = 3, include_future = TRUE)
@@ -9350,5 +9611,6 @@ usethis::use_pipe(export = TRUE)
devtools::document()
# Keep eval=FALSE to avoid infinite loop in case you hit the knit button
# Execute in the console directly
-fusen::inflate(flat_file = "dev/flat_full.Rmd", vignette_name = NA, open_vignette = F, overwrite = T)
+fusen::inflate(pkg = "/home/adminarnaud/Documents/HTRfit/", flat_file = "dev/flat_full.Rmd",
+ vignette_name = NA, open_vignette = F, overwrite = T)
```
diff --git a/man/detect_categoricals_vars.Rd b/man/detect_categoricals_vars.Rd
new file mode 100644
index 0000000..c4b7014
--- /dev/null
+++ b/man/detect_categoricals_vars.Rd
@@ -0,0 +1,32 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/update_fittedmodel.R
+\name{detect_categoricals_vars}
+\alias{detect_categoricals_vars}
+\title{Detects categorical variables based on reference labels in a glmmTMB object.}
+\usage{
+detect_categoricals_vars(tmb_frame, ref_labels)
+}
+\arguments{
+\item{tmb_frame}{The data frame of a glmmTMB object.}
+
+\item{ref_labels}{Vector of reference labels corresponding to categorical variables.}
+}
+\value{
+Names of the categorical variables detected.
+}
+\description{
+This function detects categorical variables based on reference labels in a glmmTMB object's frame.
+}
+\examples{
+data("mtcars")
+## -- specify categorical var
+mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+levels(mtcars$vs) <- c("V-shaped", "straight")
+mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+levels(mtcars$am) <- c("automatic", "manual")
+## -- For each group of number of cylinder,
+## -- explain fuel consumption with engine shape, Gross horsepower, and transmission type
+list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am,
+ data = mtcars, group_by = "cyl" , n.cores = 1)
+detect_categoricals_vars(list_tmb[["6"]]$frame, c("straight", "manual"))
+}
diff --git a/man/first_non_null_index.Rd b/man/first_non_null_index.Rd
new file mode 100644
index 0000000..1335087
--- /dev/null
+++ b/man/first_non_null_index.Rd
@@ -0,0 +1,21 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/utils.R
+\name{first_non_null_index}
+\alias{first_non_null_index}
+\title{Finds the index of the first non-null element in a list.}
+\usage{
+first_non_null_index(lst)
+}
+\arguments{
+\item{lst}{The list to search.}
+}
+\value{
+The index of the first non-null element, or NULL if no non-null element is found.
+}
+\description{
+This function searches a list and returns the index of the first non-null element.
+}
+\examples{
+my_list <- list(NULL, NULL, 3, 5, NULL)
+first_non_null_index(my_list) # Returns 3
+}
diff --git a/man/get_rsquare_2plot.Rd b/man/get_rsquare_2plot.Rd
index 8dc20c6..37bb664 100644
--- a/man/get_rsquare_2plot.Rd
+++ b/man/get_rsquare_2plot.Rd
@@ -17,6 +17,8 @@ This function takes a data frame with R-squared values,
computes position coordinates, and prepares data for plotting.
}
\examples{
-data_rsquare <- data.frame(from = c("A", "B", "C"), description = c("Desc1", "Desc2", "Desc3"), R2 = c(0.9, 0.8, 0.7))
+data_rsquare <- data.frame(from = c("A", "B", "C"),
+ description = c("Desc1", "Desc2", "Desc3"),
+ R2 = c(0.9, 0.8, 0.7))
result <- get_rsquare_2plot(data_rsquare)
}
diff --git a/man/relevel_list_tmb_frame.Rd b/man/relevel_list_tmb_frame.Rd
new file mode 100644
index 0000000..81d1c33
--- /dev/null
+++ b/man/relevel_list_tmb_frame.Rd
@@ -0,0 +1,37 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/update_fittedmodel.R
+\name{relevel_list_tmb_frame}
+\alias{relevel_list_tmb_frame}
+\title{Re-levels categorical variables in a the frame of a list of glmmTMB objects.}
+\usage{
+relevel_list_tmb_frame(list_tmb, categorical_vars, ref_labels)
+}
+\arguments{
+\item{list_tmb}{List of glmmTMB objects.}
+
+\item{categorical_vars}{Names of the categorical variables to be re-leveled.}
+
+\item{ref_labels}{Vector of reference labels corresponding to the categorical variables.}
+}
+\value{
+A list of glmmTMB objects with re-leveled categorical variables.
+}
+\description{
+This function re-levels categorical variables in a list of glmmTMB objects using the specified reference labels.
+}
+\examples{
+# Example usage:
+# -- Load the mtcars dataset
+data("mtcars")
+## -- specify categorical var
+mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+levels(mtcars$vs) <- c("V-shaped", "straight")
+mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+levels(mtcars$am) <- c("automatic", "manual")
+# -- For each group of number of cylinders,
+# -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am,
+ data = mtcars, group_by = "cyl", n.cores = 1)
+# -- Relevel transmission and engine shape variables
+relevel_list_tmb_frame(list_tmb, c("am", "vs"), c("manual", "straight"))
+}
diff --git a/man/relevelling_factors.Rd b/man/relevelling_factors.Rd
new file mode 100644
index 0000000..94013d0
--- /dev/null
+++ b/man/relevelling_factors.Rd
@@ -0,0 +1,32 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/update_fittedmodel.R
+\name{relevelling_factors}
+\alias{relevelling_factors}
+\title{Relevels factors in a list of glmmTMB objects using specified reference labels.}
+\usage{
+relevelling_factors(list_tmb, reference_labels)
+}
+\arguments{
+\item{list_tmb}{List of glmmTMB objects.}
+
+\item{reference_labels}{Vector of reference labels.}
+}
+\value{
+A list of glmmTMB objects with re-leveled factors.
+}
+\description{
+This function re-levels factors in a list of glmmTMB objects based on the specified reference labels.
+}
+\examples{
+data("mtcars")
+## -- specify categorical var
+mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+levels(mtcars$vs) <- c("V-shaped", "straight")
+mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+levels(mtcars$am) <- c("automatic", "manual")
+## -- For each group of number of cylinder,
+## -- explain fuel consumption with engine shape, Gross horsepower, and transmission type
+list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am,
+ data = mtcars, group_by = "cyl", n.cores = 1 )
+relevelling_factors(list_tmb , c("straight", "manual"))
+}
diff --git a/man/updateParallel.Rd b/man/updateParallel.Rd
index 665857d..616c6cc 100644
--- a/man/updateParallel.Rd
+++ b/man/updateParallel.Rd
@@ -7,6 +7,7 @@
updateParallel(
formula,
list_tmb,
+ reference_labels = c(),
n.cores = NULL,
cl_type = "PSOCK",
log_file = paste(tempdir(check = FALSE), "htrfit.log", sep = "/"),
@@ -18,9 +19,11 @@ updateParallel(
\item{list_tmb}{List of glmmTMB objects.}
+\item{reference_labels}{Vector of reference labels. Default is c(), selecting the first alphanumeric label as reference.}
+
\item{n.cores}{Number of cores to use for parallel processing. If NULL, the function will use all available cores.}
-\item{cl_type}{cluster type (defautl "PSOCK"). "FORK" is recommanded for linux.}
+\item{cl_type}{cluster type (default "PSOCK"). "FORK" is recommended for linux.}
\item{log_file}{File path for the log output (default: Rtmpdir/htrfit.log).}
@@ -30,9 +33,12 @@ updateParallel(
A list of updated GLMNB models.
}
\description{
-This function fits glmmTMB models in parallel using multiple cores, allowing for faster computation.
+This function updates glmmTMB models in parallel using multiple cores, allowing for faster computation.
+It updates the models with new reference labels if specified.
+It can also be used to fit a new formula or to change additional parameters of glmmTMB (param : "...").
}
\examples{
+# -- Example usage: update formula
data(iris)
groups <- unique(iris$Species)
group_by <- "Species"
@@ -40,4 +46,19 @@ formula <- Sepal.Length ~ Sepal.Width + Petal.Length
fitted_models <- fitModelParallel(formula, iris, group_by, n.cores = 1)
new_formula <- Sepal.Length ~ Sepal.Width
results <- updateParallel(new_formula, fitted_models, n.cores = 1)
+# Example usage: update reference
+# -- Load the mtcars dataset
+data("mtcars")
+# -- Specify categorical variables
+mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+levels(mtcars$vs) <- c("V-shaped", "straight")
+mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+levels(mtcars$am) <- c("automatic", "manual")
+# -- For each group of number of cylinders:
+# -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am,
+ data = mtcars, group_by = "cyl", n.cores = 1)
+# -- Relevel transmission and engine shape variables
+list_tmb <- updateParallel(formula = mpg ~ hp + vs + am, list_tmb,
+ reference_labels = c("straight", "manual"), n.cores = 1)
}
diff --git a/tests/testthat/test-actual_interactionfixeffects.R b/tests/testthat/test-actual_interactionfixeffects.R
index 3dea605..fdd12c2 100644
--- a/tests/testthat/test-actual_interactionfixeffects.R
+++ b/tests/testthat/test-actual_interactionfixeffects.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
test_that("filter_dataframe retourne le dataframe filtré correctement", {
diff --git a/tests/testthat/test-actual_mainfixeffects.R b/tests/testthat/test-actual_mainfixeffects.R
index 0c47b6a..308ad81 100644
--- a/tests/testthat/test-actual_mainfixeffects.R
+++ b/tests/testthat/test-actual_mainfixeffects.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
test_that("Test for subsetFixEffectInferred function", {
diff --git a/tests/testthat/test-anova.R b/tests/testthat/test-anova.R
index 5959558..d8b544b 100644
--- a/tests/testthat/test-anova.R
+++ b/tests/testthat/test-anova.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/tests/testthat/test-basal_expression_scaling.R b/tests/testthat/test-basal_expression_scaling.R
index 86876f8..170623f 100644
--- a/tests/testthat/test-basal_expression_scaling.R
+++ b/tests/testthat/test-basal_expression_scaling.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
test_that("generate_basal_expression returns correct number of genes", {
diff --git a/tests/testthat/test-counts_plot.R b/tests/testthat/test-counts_plot.R
index e0c29ee..5f2bad2 100644
--- a/tests/testthat/test-counts_plot.R
+++ b/tests/testthat/test-counts_plot.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/tests/testthat/test-datafrommvrnorm_manipulations.R b/tests/testthat/test-datafrommvrnorm_manipulations.R
index 7e1bb39..8b2dcbb 100644
--- a/tests/testthat/test-datafrommvrnorm_manipulations.R
+++ b/tests/testthat/test-datafrommvrnorm_manipulations.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
test_that("getInput2mvrnorm returns the correct list", {
list_var <- init_variable()
diff --git a/tests/testthat/test-datafromuser_manipulations.R b/tests/testthat/test-datafromuser_manipulations.R
index f2fd8db..a8997d1 100644
--- a/tests/testthat/test-datafromuser_manipulations.R
+++ b/tests/testthat/test-datafromuser_manipulations.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
# Test unitaires pour la fonction join_dtf
test_that("join_dtf réalise la jointure correctement", {
diff --git a/tests/testthat/test-evaluate_dispersion.R b/tests/testthat/test-evaluate_dispersion.R
index 46c4bce..cb72158 100644
--- a/tests/testthat/test-evaluate_dispersion.R
+++ b/tests/testthat/test-evaluate_dispersion.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/tests/testthat/test-evaluation_identity.R b/tests/testthat/test-evaluation_identity.R
index 1dedef1..793bee5 100644
--- a/tests/testthat/test-evaluation_identity.R
+++ b/tests/testthat/test-evaluation_identity.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/tests/testthat/test-evaluation_withmixedeffect.R b/tests/testthat/test-evaluation_withmixedeffect.R
index fbaee4d..0d7af6f 100644
--- a/tests/testthat/test-evaluation_withmixedeffect.R
+++ b/tests/testthat/test-evaluation_withmixedeffect.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/tests/testthat/test-fitmodel.R b/tests/testthat/test-fitmodel.R
index e6a11fc..5b0d2e8 100644
--- a/tests/testthat/test-fitmodel.R
+++ b/tests/testthat/test-fitmodel.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
@@ -90,8 +90,9 @@ test_that("Detect rank-deficient model matrix and throw error", {
z = factor(rep(c("zA","zB"), each = 5)),
y = rnorm(10))
formula <- y ~ x + w + z + y:w
- expect_error(is_fullrank(metadata, formula),
- regexp = "The model matrix is not full rank, so the model cannot be fit as specified.")
+ expect_warning(is_fullrank(metadata, formula))
+ res <- suppressWarnings(is_fullrank(metadata, formula))
+ expect_false(res)
})
# Test if a rank-deficient model matrix is detected and throws an error
@@ -101,8 +102,9 @@ test_that("Detect rank-deficient model matrix and throw error (with random eff)"
z = factor(rep(c("zA","zB"), each = 5)),
y = rnorm(10))
formula <- y ~ x + w + z + y:w + (1 | w)
- expect_error(is_fullrank(metadata, formula),
- regexp = "The model matrix is not full rank, so the model cannot be fit as specified.")
+ expect_warning(is_fullrank(metadata, formula))
+ res <- suppressWarnings(is_fullrank(metadata, formula))
+ expect_false(res)
})
# Test if a rank-deficient model matrix is detected and throws an error
diff --git a/tests/testthat/test-glance_glmmtmb.R b/tests/testthat/test-glance_glmmtmb.R
index 9e5d6cc..6b1c4a4 100644
--- a/tests/testthat/test-glance_glmmtmb.R
+++ b/tests/testthat/test-glance_glmmtmb.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
test_that("glance_tmb returns the summary statistics for multiple models", {
diff --git a/tests/testthat/test-mock_rnaseq.R b/tests/testthat/test-mock_rnaseq.R
index ac3bcf3..6596bc4 100644
--- a/tests/testthat/test-mock_rnaseq.R
+++ b/tests/testthat/test-mock_rnaseq.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
# Test for is_dispersionMatrixValid
diff --git a/tests/testthat/test-plot_metrics.R b/tests/testthat/test-plot_metrics.R
index d52f49d..777d383 100644
--- a/tests/testthat/test-plot_metrics.R
+++ b/tests/testthat/test-plot_metrics.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/tests/testthat/test-precision_recall.R b/tests/testthat/test-precision_recall.R
index da44700..047b1b7 100644
--- a/tests/testthat/test-precision_recall.R
+++ b/tests/testthat/test-precision_recall.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/tests/testthat/test-prepare_data2fit.R b/tests/testthat/test-prepare_data2fit.R
index beb3720..4e7a7b6 100644
--- a/tests/testthat/test-prepare_data2fit.R
+++ b/tests/testthat/test-prepare_data2fit.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/tests/testthat/test-receiver_operating_characteristic.R b/tests/testthat/test-receiver_operating_characteristic.R
index 4e4aeac..016cf0a 100644
--- a/tests/testthat/test-receiver_operating_characteristic.R
+++ b/tests/testthat/test-receiver_operating_characteristic.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/tests/testthat/test-sequencing_depth_scaling.R b/tests/testthat/test-sequencing_depth_scaling.R
index 0f1e0a8..294f05f 100644
--- a/tests/testthat/test-sequencing_depth_scaling.R
+++ b/tests/testthat/test-sequencing_depth_scaling.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
# Test case 1: Scaling with valid min_seq_depth and max_seq_depth
diff --git a/tests/testthat/test-setcorrelation.R b/tests/testthat/test-setcorrelation.R
index e70b19d..118479d 100644
--- a/tests/testthat/test-setcorrelation.R
+++ b/tests/testthat/test-setcorrelation.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
test_that("compute_covariation returns the correct covariation", {
diff --git a/tests/testthat/test-simulation.R b/tests/testthat/test-simulation.R
index c39d27f..0868462 100644
--- a/tests/testthat/test-simulation.R
+++ b/tests/testthat/test-simulation.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
# Test case 1: Check if the function returns a data frame
test_that("getInput2simulation returns a data frame", {
diff --git a/tests/testthat/test-simulation_initialization.R b/tests/testthat/test-simulation_initialization.R
index 01dffa2..6c3263b 100644
--- a/tests/testthat/test-simulation_initialization.R
+++ b/tests/testthat/test-simulation_initialization.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
test_that("endsWithDigit returns the correct result", {
diff --git a/tests/testthat/test-simulation_report.R b/tests/testthat/test-simulation_report.R
index 4f06316..cdca5d2 100644
--- a/tests/testthat/test-simulation_report.R
+++ b/tests/testthat/test-simulation_report.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
diff --git a/tests/testthat/test-subsetgenes.R b/tests/testthat/test-subsetgenes.R
index 40ebcb8..e79d2fe 100644
--- a/tests/testthat/test-subsetgenes.R
+++ b/tests/testthat/test-subsetgenes.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
test_that("subsetGenes return correct ouptut", {
N_GENES = 100
diff --git a/tests/testthat/test-tidy_glmmtmb.R b/tests/testthat/test-tidy_glmmtmb.R
index 06ce9c4..40aff2d 100644
--- a/tests/testthat/test-tidy_glmmtmb.R
+++ b/tests/testthat/test-tidy_glmmtmb.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
test_that("extract_fixed_effect returns the correct results for glmmTMB models", {
diff --git a/tests/testthat/test-update_fittedmodel.R b/tests/testthat/test-update_fittedmodel.R
index 74702cd..9bdadb0 100644
--- a/tests/testthat/test-update_fittedmodel.R
+++ b/tests/testthat/test-update_fittedmodel.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
# Test updateParallel function
test_that("updateParallel function returns correct results", {
@@ -42,8 +42,105 @@ test_that("updateParallel function returns correct results", {
n.cores = 1,
family = glmmTMB::ziGamma(link = "inverse")))
expect_s3_class(updated_updated_model$setosa$call$family, "family")
+
+
+ ## -- update label reference
+ data("mtcars")
+ # -- Specify categorical variables
+ mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+ levels(mtcars$vs) <- c("V-shaped", "straight")
+ mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+ levels(mtcars$am) <- c("automatic", "manual")
+ # -- For each group of number of cylinders:
+ # -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+ list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am, data = mtcars, group_by = "cyl", n.cores = 1)
+ # -- Relevel transmission and engine shape variables
+ reference_labels <- c("straight", "manual")
+ result <- updateParallel(formula = mpg ~ hp + vs + am, list_tmb = list_tmb, reference_labels = reference_labels, n.cores = 1)
+
+ # Check if the returned list has the same length as the mock list
+ expect_equal(length(result), length(list_tmb))
+ expect_equal(levels(result[["6"]]$frame$am)[1] , "manual")
+ expect_equal(levels(result[["4"]]$frame$am)[1] , "manual")
+ expect_equal(levels(result[["6"]]$frame$vs)[1] , "straight")
+ expect_equal(levels(result[["4"]]$frame$vs)[1] , "straight")
+
+})
+
+
+
+# Test for relevel_list_tmb_frame function
+test_that("relevel_list_tmb_frame re-levels categorical variables correctly", {
+ data("mtcars")
+ # -- Specify categorical variables
+ mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+ levels(mtcars$vs) <- c("V-shaped", "straight")
+ mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+ levels(mtcars$am) <- c("automatic", "manual")
+ # -- For each group of number of cylinders:
+ # -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+ list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am, data = mtcars, group_by = "cyl", n.cores = 1)
+ # -- Relevel transmission and engine shape variables
+ reference_labels <- c("straight", "manual")
+ result <- relevel_list_tmb_frame(list_tmb, c("vs", "am"), reference_labels)
+
+ # Check if the returned list has the same length as the mock list
+ expect_equal(length(result), length(list_tmb))
+ # Check if categorical variables have been re-leveled correctly
+ expect_equal(levels(result[["6"]]$frame$am)[1] , "manual")
+ expect_equal(levels(result[["4"]]$frame$am)[1] , "manual")
+ expect_equal(levels(result[["6"]]$frame$vs)[1] , "straight")
+ expect_equal(levels(result[["4"]]$frame$vs)[1] , "straight")
})
+# Test for relevelling_factors function
+test_that("relevelling_factors relevel categorical variables correctly", {
+ data("mtcars")
+ # -- Specify categorical variables
+ mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+ levels(mtcars$vs) <- c("V-shaped", "straight")
+ mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+ levels(mtcars$am) <- c("automatic", "manual")
+ # -- For each group of number of cylinders:
+ # -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+ list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am, data = mtcars, group_by = "cyl", n.cores = 1)
+ # -- Relevel transmission and engine shape variables
+ reference_labels <- c("straight", "manual")
+ result <- relevelling_factors(list_tmb, reference_labels)
+
+ # Check if the returned list has the same length as the mock list
+ expect_equal(length(result), length(list_tmb))
+ # Check if categorical variables have been re-leveled correctly
+ expect_equal(levels(result[["6"]]$frame$am)[1] , "manual")
+ expect_equal(levels(result[["4"]]$frame$am)[1] , "manual")
+ expect_equal(levels(result[["6"]]$frame$vs)[1] , "straight")
+ expect_equal(levels(result[["4"]]$frame$vs)[1] , "straight")
+})
+
+
+
+# Test for detect_categoricals_vars function
+test_that("detect_categoricals_vars detect categorical variables correctly", {
+ data("mtcars")
+ # -- Specify categorical variables
+ mtcars$vs <- factor(mtcars$vs) ## Engine (0 = V-shaped, 1 = straight)
+ levels(mtcars$vs) <- c("V-shaped", "straight")
+ mtcars$am <- factor(mtcars$am) ## Transmission (0 = automatic, 1 = manual)
+ levels(mtcars$am) <- c("automatic", "manual")
+ # -- For each group of number of cylinders:
+ # -- Explain fuel consumption with engine shape, Gross horsepower, and transmission type
+ list_tmb <- fitModelParallel(formula = mpg ~ hp + vs + am, data = mtcars, group_by = "cyl", n.cores = 1)
+ # -- Relevel transmission and engine shape variables
+ reference_labels <- c("straight", "manual")
+ result <- detect_categoricals_vars(list_tmb[["6"]]$frame, reference_labels)
+
+ # Check if the returned list has the same length as the mock list
+ expect_equal(length(result), length(reference_labels))
+ # Check if categorical variables have been re-leveled correctly
+ expect_equal(result, c("vs", "am"))
+})
+
+
# Test parallel_update function
test_that("parallel_update function returns correct results", {
# Load the required data
diff --git a/tests/testthat/test-utils.R b/tests/testthat/test-utils.R
index bca87e2..dbf6fa0 100644
--- a/tests/testthat/test-utils.R
+++ b/tests/testthat/test-utils.R
@@ -1,4 +1,13 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
+
+
+# Test for first_non_null_index function
+test_that("first_non_null_index returns the correct index", {
+ lst <- list(NULL, NULL, 3, 5, NULL)
+ expect_equal(first_non_null_index(lst), 3)
+})
+
+
# Test unitaires pour la fonction join_dtf
test_that("join_dtf réalise la jointure correctement", {
diff --git a/tests/testthat/test-waldtest.R b/tests/testthat/test-waldtest.R
index 9aad7da..75d872e 100644
--- a/tests/testthat/test-waldtest.R
+++ b/tests/testthat/test-waldtest.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
# Test unitaires
diff --git a/tests/testthat/test-wrapper_dds.R b/tests/testthat/test-wrapper_dds.R
index 03d57e1..59a2d4e 100644
--- a/tests/testthat/test-wrapper_dds.R
+++ b/tests/testthat/test-wrapper_dds.R
@@ -1,4 +1,4 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
+# WARNING - Generated by {fusen} from dev/flat_full.Rmd: do not edit by hand
@@ -119,3 +119,4 @@ test_that("wrapperDESeq2 function works correctly", {
})
+
diff --git a/vignettes/03-rnaseq_analysis.Rmd b/vignettes/03-rnaseq_analysis.Rmd
index bbbc423..e87f7fd 100644
--- a/vignettes/03-rnaseq_analysis.Rmd
+++ b/vignettes/03-rnaseq_analysis.Rmd
@@ -139,7 +139,7 @@ my_tidy_res[1:3,]
## Update fit
-The `updateParallel()` function updates and re-fits a model for each gene. It offers options similar to those in `fitModelParallel()`.
+The `updateParallel()` function updates and re-fits a model for each gene. It offers options similar to those in `fitModelParallel()`. In addition, it is possible to modify the reference level of the categorical variable used in your model in order to use different contrast.
```{r example-update, warning = FALSE, message = FALSE}
@@ -166,6 +166,16 @@ l_tmb <- updateParallel(
list_tmb = l_tmb ,
family = glmmTMB::nbinom2(link = "log"),
n.cores = 1)
+
+## -- modif reference levels
+## -- genotype reference = "genotype2"
+## -- environment reference = "environment2"
+l_tmb <- updateParallel(
+ formula = kij ~ genotype + environment ,
+ list_tmb = l_tmb ,
+ family = glmmTMB::nbinom2(link = "log"),
+ n.cores = 1,
+ reference_labels = c("genotype2", "environment2"))
```
#### Struture of list tmb object
--
GitLab