diff --git a/NAMESPACE b/NAMESPACE
index af9a2dc9c5b15b6b7dc17debdfede65b7fff193f..97817f4ab418db8446ef1b0909d5ba6011ab3f4d 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -1,15 +1,6 @@
 # Generated by roxygen2: do not edit by hand
 
 export("%>%")
-export(.fitModel)
-export(.getColumnWithSampleID)
-export(.isDispersionMatrixValid)
-export(.parallel_fit)
-export(.parallel_update)
-export(.replicateByGroup)
-export(.replicateMatrix)
-export(.replicateRows)
-export(.subsetData_andfit)
 export(addBasalExpression)
 export(add_interaction)
 export(already_init_variable)
@@ -33,15 +24,16 @@ export(drop_randfx)
 export(endsWithDigit)
 export(evaluateDispersion)
 export(exportReportFile)
-export(extractDESeqDispersion)
-export(extractTMBDispersion)
+export(extract_ddsDispersion)
 export(extract_fixed_effect)
 export(extract_ran_pars)
+export(extract_tmbDispersion)
 export(fillInCovarMatrice)
 export(fillInInteraction)
 export(fillInVariable)
 export(filter_dataframe)
 export(findAttribute)
+export(fitModel)
 export(fitModelParallel)
 export(fitUpdate)
 export(generateActualForMainFixEff)
@@ -50,7 +42,7 @@ export(generateActualInteractionX3_FixEff)
 export(generateCountTable)
 export(generateGridCombination_fromListVar)
 export(generateReplicationMatrix)
-export(generate_BE)
+export(generate_basal_expression)
 export(getActualInteractionFixEff)
 export(getActualIntercept)
 export(getActualMainFixEff)
@@ -58,6 +50,7 @@ export(getActualMixed_typeI)
 export(getBinExpression)
 export(getCategoricalVar_inFixedEffect)
 export(getCoefficients)
+export(getColumnWithSampleID)
 export(getCovarianceMatrix)
 export(getData2computeActualFixEffect)
 export(getDataFromMvrnorm)
@@ -66,12 +59,15 @@ export(getDispersionComparison)
 export(getDispersionMatrix)
 export(getEstimate_df)
 export(getGeneMetadata)
+export(getGivenAttribute)
 export(getGlance)
 export(getGridCombination)
 export(getGrobTable)
 export(getInput2mvrnorm)
 export(getInput2simulation)
 export(getLabelExpected)
+export(getLabels)
+export(getListVar)
 export(getLog_qij)
 export(getMu_ij)
 export(getMu_ij_matrix)
@@ -84,7 +80,7 @@ export(getSettingsTable)
 export(getStandardDeviationInCorrelation)
 export(getTidyGlmmTMB)
 export(getValidDispersion)
-export(get_inference)
+export(get_inference_dds)
 export(glance_tmb)
 export(group_logQij_per_genes_and_labels)
 export(handleAnovaError)
@@ -94,6 +90,7 @@ export(inferenceToExpected_withMixedEff)
 export(init_variable)
 export(inputs_checking)
 export(isValidInput2fit)
+export(is_dispersionMatrixValid)
 export(is_formula_mixedTypeI)
 export(is_fullrank)
 export(is_mixedEffect_inFormula)
@@ -104,18 +101,24 @@ export(launchUpdate)
 export(medianRatioNormalization)
 export(metrics_plot)
 export(mock_rnaseq)
+export(parallel_fit)
+export(parallel_update)
 export(prepareData2computeInteraction)
 export(prepareData2fit)
 export(removeDigitsAtEnd)
 export(removeDuplicatedWord)
 export(renameColumns)
 export(reorderColumns)
+export(replicateByGroup)
+export(replicateMatrix)
+export(replicateRows)
 export(roc_plot)
 export(samplingFromMvrnorm)
 export(scaleCountsTable)
 export(set_correlation)
 export(simulationReport)
 export(subsetByTermLabel)
+export(subsetData_andfit)
 export(subsetFixEffectInferred)
 export(subsetGenes)
 export(subset_glance)
@@ -123,8 +126,9 @@ export(tidy_results)
 export(tidy_tmb)
 export(updateParallel)
 export(wald_test)
-export(wrapper_DESeq2)
+export(wrap_dds)
 export(wrapper_var_cor)
+importFrom(MASS,mvrnorm)
 importFrom(car,Anova)
 importFrom(data.table,data.table)
 importFrom(data.table,setDT)
@@ -161,6 +165,7 @@ importFrom(rlang,":=")
 importFrom(stats,anova)
 importFrom(stats,as.formula)
 importFrom(stats,cor)
+importFrom(stats,drop.terms)
 importFrom(stats,median)
 importFrom(stats,model.matrix)
 importFrom(stats,p.adjust)
diff --git a/R/actualinteractionfixeffects.R b/R/actual_interactionfixeffects.R
similarity index 100%
rename from R/actualinteractionfixeffects.R
rename to R/actual_interactionfixeffects.R
diff --git a/R/actualmainfixeffects.R b/R/actual_mainfixeffects.R
similarity index 94%
rename from R/actualmainfixeffects.R
rename to R/actual_mainfixeffects.R
index 09ecd4a593403f7d5d067a8654ba1b99d771dca3..6642c14a199445e3c2a2c6dd7bb1fcf06b195c46 100644
--- a/R/actualmainfixeffects.R
+++ b/R/actual_mainfixeffects.R
@@ -19,22 +19,6 @@ averageByGroup <- function(data, column, group_by) {
   return(result)
 }
 
-#' Convert specified columns to factor
-#'
-#' @param data The input data frame
-#' @param columns The column names to be converted to factors
-#' @return The modified data frame with specified columns converted to factors
-#' @export
-convert2Factor <- function(data, columns) {
-  if (is.character(columns)) {
-    columns <- match(columns, colnames(data))
-  }
-
-  if (length(columns) > 1) data[, columns] <- lapply(data[, columns], as.factor )
-  else data[, columns] <- as.factor(data[, columns])
-  return(data)
-}
-
 #' Subset Fixed Effect Inferred Terms
 #'
 #' This function subsets the tidy TMB object to extract the fixed effect inferred terms
@@ -240,7 +224,3 @@ getActualMainFixEff <- function( l_term , fixeEff_dataActual , df_actualIntercep
   return(df_actual)
 }
 
-
-
-
-
diff --git a/R/anova.R b/R/anova.R
index ce7c023ae3e36331d1fe7e72647183656d0056fe..c2a2c75b93412a58981f6bd3b3b2c2e050e43ffa 100644
--- a/R/anova.R
+++ b/R/anova.R
@@ -5,7 +5,7 @@
 #'
 #' This function handles ANOVA errors and warnings during the ANOVA calculation process.
 #'
-#' @param l_TMB A list of fitted glmmTMB models.
+#' @param list_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.
 #' 
@@ -13,17 +13,17 @@
 #' @export
 #' 
 #' @examples
-#' l_tmb <- fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length,
+#' list_tmb <- fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length,
 #'                           data = iris, group_by = "Species", n.cores = 1)
-#' anova_res <- handleAnovaError(l_tmb, "setosa", type = "III")
+#' anova_res <- handleAnovaError(list_tmb, "setosa", type = "III")
 #'
 #' @importFrom car Anova
 #' @export
-handleAnovaError <- function(l_TMB, group, ...) {
+handleAnovaError <- function(list_tmb, group, ...) {
   tryCatch(
     expr = {
       withCallingHandlers(
-        car::Anova(l_TMB[[group]], ...),
+        car::Anova(list_tmb[[group]], ...),
         warning = function(w) {
           message(paste(Sys.time(), "warning for group", group, ":", conditionMessage(w)))
           invokeRestart("muffleWarning")
@@ -43,7 +43,7 @@ handleAnovaError <- function(l_TMB, group, ...) {
 #' models in parallel for different groups specified in the list. It returns a list
 #' of ANOVA results for each group.
 #'
-#' @param l_tmb A list of \code{glmmTMB} models, with model names corresponding to the groups.
+#' @param list_tmb 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.
@@ -51,14 +51,14 @@ handleAnovaError <- function(l_TMB, group, ...) {
 #' @examples
 #' # Perform ANOVA
 #' data(iris)
-#' l_tmb<- fitModelParallel( Sepal.Length ~ Sepal.Width  + Petal.Length, 
+#' list_tmb<- fitModelParallel( Sepal.Length ~ Sepal.Width  + Petal.Length, 
 #'                          data = iris, group_by = "Species", n.cores = 1 )
-#' anov_res <- anovaParallel(l_tmb , type = "III")
+#' anov_res <- anovaParallel(list_tmb , type = "III")
 #' @importFrom stats anova
 #' @export
-anovaParallel <- function(l_tmb, ...) {
-  l_group <- attributes(l_tmb)$names
-  lapply(stats::setNames(l_group, l_group), function(group) handleAnovaError(l_tmb, group, ...))
+anovaParallel <- function(list_tmb, ...) {
+  l_group <- attributes(list_tmb)$names
+  lapply(stats::setNames(l_group, l_group), function(group) handleAnovaError(list_tmb, group, ...))
 }
 
 
diff --git a/R/scalinggeneexpression.R b/R/basal_expression__scaling.R
similarity index 87%
rename from R/scalinggeneexpression.R
rename to R/basal_expression__scaling.R
index 19729b5c82a831e451830b8dd6630044ce0ce749..c2ddc5d3bbcaf8167cc73024562309602443501d 100644
--- a/R/scalinggeneexpression.R
+++ b/R/basal_expression__scaling.R
@@ -30,7 +30,7 @@ getBinExpression <- function(dtf_coef, n_bins){
 
 #' Generate BE data.
 #' 
-#' This function generates BE data for a given number of genes, in a vector of BE values.
+#' This function generates basal expression data for a given number of genes, in a vector of basal expression values.
 #' 
 #' @param n_genes The number of genes to generate BE data for.
 #' @param basal_expression a numeric vector from which sample BE for eacg genes
@@ -38,10 +38,10 @@ getBinExpression <- function(dtf_coef, n_bins){
 #' @return A data frame containing gene IDs, BE values
 #' 
 #' @examples
-#' generate_BE(n_genes = 100, 10)
+#' generate_basal_expression(n_genes = 100, 10)
 #' 
 #' @export
-generate_BE <- function(n_genes, basal_expression) {
+generate_basal_expression <- function(n_genes, basal_expression) {
   ## --avoid bug if one value in basal_expr
   pool2sample <- c(basal_expression, basal_expression)
   BE <- sample(x = pool2sample, size = n_genes, replace = T)
@@ -56,7 +56,7 @@ generate_BE <- function(n_genes, basal_expression) {
 #'
 #' This function takes the coefficients data frame \code{dtf_coef} and computes
 #' basal expression for gene expression. The scaling factors are generated 
-#' using the function \code{generate_BE}.
+#' using the function \code{generate_basal_expression}.
 #'
 #' @param dtf_coef A data frame containing the coefficients for gene expression.
 #' @param n_genes number of genes in simulation
@@ -72,7 +72,7 @@ generate_BE <- function(n_genes, basal_expression) {
 #' dtf_coef <- getLog_qij(dtf_coef)
 #' addBasalExpression(dtf_coef, N_GENES, 1)
 addBasalExpression <- function(dtf_coef, n_genes, basal_expression){
-    BE_df  <-  generate_BE(n_genes, basal_expression )
+    BE_df  <-  generate_basal_expression(n_genes, basal_expression )
     dtf_coef <- join_dtf(dtf_coef, BE_df, "geneID", "geneID")
     return(dtf_coef) 
 }
diff --git a/R/countsplot.R b/R/counts_plot.R
similarity index 100%
rename from R/countsplot.R
rename to R/counts_plot.R
diff --git a/R/datafrommvrnorm_manipulations.R b/R/datafrommvrnorm_manipulations.R
index 5b9ae564f0c33c43f7380f911e86d7dd03df570a..cbc046f574486dede578eadb8adce36e329815d7 100644
--- a/R/datafrommvrnorm_manipulations.R
+++ b/R/datafrommvrnorm_manipulations.R
@@ -94,37 +94,6 @@ fillInCovarMatrice <- function(covarMatrice, covar){
   return(covarMatrice)
 }
 
-
-#' Check if a matrix is positive definite
-#' This function checks whether a given matrix is positive definite, i.e., all of its eigenvalues are positive.
-#' @param mat The matrix to be checked.
-#' @return A logical value indicating whether the matrix is positive definite.
-#' @export
-#' @examples
-#' # Create a positive definite matrix
-#' mat1 <- matrix(c(4, 2, 2, 3), nrow = 2)
-#' is_positive_definite(mat1)
-#' # Expected output: TRUE
-#'
-#' # Create a non-positive definite matrix
-#' mat2 <- matrix(c(4, 2, 2, -3), nrow = 2)
-#' is_positive_definite(mat2)
-#' # Expected output: FALSE
-#'
-#' # Check an empty matrix
-#' mat3 <- matrix(nrow = 0, ncol = 0)
-#' is_positive_definite(mat3)
-#' # Expected output: TRUE
-#'
-#' @export
-is_positive_definite <- function(mat) {
-  if (nrow(mat) == 0 && ncol(mat) == 0) return(TRUE)
-  eigenvalues <- eigen(mat)$values
-  all(eigenvalues > 0)
-}
-
-
-
 #' getGeneMetadata
 #'
 #' @inheritParams init_variable
@@ -141,7 +110,7 @@ is_positive_definite <- function(mat) {
 getGeneMetadata <- function(list_var, n_genes) {
   metaData <- generateGridCombination_fromListVar(list_var)
   n_combinations <- dim(metaData)[1]
-  genes_vec <- base::paste("gene", 1:n_genes, sep = "")
+  genes_vec <- paste("gene", 1:n_genes, sep = "")
   geneID <- rep(genes_vec, each = n_combinations)
   metaData <- cbind(geneID, metaData)
   
@@ -190,7 +159,7 @@ getDataFromMvrnorm <- function(list_var, input2mvrnorm, n_genes = 1) {
 #' samples generated from multivariate normal distribution
 #' 
 #' @export
-#'
+#' @importFrom MASS mvrnorm
 #' @examples
 #' n <- 100
 #' mu <- c(0, 0)
@@ -198,7 +167,6 @@ getDataFromMvrnorm <- function(list_var, input2mvrnorm, n_genes = 1) {
 #' samples <- samplingFromMvrnorm(n_samplings = n, l_mu = mu, matx_cov = covMatrix)
 samplingFromMvrnorm <- function(n_samplings, l_mu, matx_cov) {
   mvrnormSamp <-  MASS::mvrnorm(n = n_samplings, mu = l_mu, Sigma = matx_cov, empirical = TRUE)
-  
   return(mvrnormSamp)
 }
 
diff --git a/R/evaluatedispersion.R b/R/evaluate_dispersion.R
similarity index 87%
rename from R/evaluatedispersion.R
rename to R/evaluate_dispersion.R
index 39b81b9f877fbd1fdd43775a684d84a78d8b7c0e..a16f96a514fd8ae6760bc4fafadf10627fae5cc2 100644
--- a/R/evaluatedispersion.R
+++ b/R/evaluate_dispersion.R
@@ -40,7 +40,6 @@ evaluateDispersion <- function(TMB_dispersion_df, DESEQ_dispersion_df, color2use
 #' @examples
 #' \dontrun{
 #' dispersion_comparison <- getDispersionComparison(inferred_disp, actual_disp)
-#' print(dispersion_comparison)
 #' }
 getDispersionComparison <- function(inferred_dispersion, actual_dispersion) {
   actual_disp <- data.frame(actual_dispersion = actual_dispersion)
@@ -55,7 +54,7 @@ getDispersionComparison <- function(inferred_dispersion, actual_dispersion) {
 #'
 #' Extracts inferred dispersion values from a DESeq2 wrapped object.
 #'
-#' @param deseq_wrapped A DESeq2 wrapped object containing dispersion values.
+#' @param dds_wrapped A DESeq2 wrapped object containing dispersion values.
 #'
 #' @return A data frame containing inferred dispersion values.
 #' 
@@ -63,11 +62,10 @@ getDispersionComparison <- function(inferred_dispersion, actual_dispersion) {
 #'
 #' @examples
 #' \dontrun{
-#' dispersion_df <- extractDESeqDispersion(deseq2_object)
-#' print(dispersion_df)
+#' dispersion_df <- extract_ddsDispersion(deseq2_object)
 #' }
-extractDESeqDispersion <- function(deseq_wrapped) {
-  inferred_dispersion <- data.frame(inferred_dispersion = deseq_wrapped$dispersion)
+extract_ddsDispersion <- function(dds_wrapped) {
+  inferred_dispersion <- data.frame(inferred_dispersion = dds_wrapped$dispersion)
   inferred_dispersion$geneID <- rownames(inferred_dispersion)
   rownames(inferred_dispersion) <- NULL
   return(inferred_dispersion)
@@ -78,7 +76,7 @@ extractDESeqDispersion <- function(deseq_wrapped) {
 #'
 #' Extracts inferred dispersion values from a TMB result object.
 #'
-#' @param l_tmb A TMB result object containing dispersion values.
+#' @param list_tmb A TMB result object containing dispersion values.
 #'
 #' @return A data frame containing inferred dispersion values.
 #' 
@@ -86,11 +84,10 @@ extractDESeqDispersion <- function(deseq_wrapped) {
 #'
 #' @examples
 #' \dontrun{
-#' dispersion_df <- extractTMBDispersion(tmb_result)
-#' print(dispersion_df)
+#' dispersion_df <- extract_tmbDispersion(tmb_result)
 #' }
-extractTMBDispersion <- function(l_tmb) {
-  glanceRes <- glance_tmb(l_tmb)
+extract_tmbDispersion <- function(list_tmb) {
+  glanceRes <- glance_tmb(list_tmb)
   inferred_dispersion <- data.frame(inferred_dispersion = glanceRes$dispersion)
   inferred_dispersion$geneID <- rownames(glanceRes)
   rownames(inferred_dispersion) <- NULL
diff --git a/R/evaluationwithmixedeffect.R b/R/evaluation_withmixedeffect.R
similarity index 100%
rename from R/evaluationwithmixedeffect.R
rename to R/evaluation_withmixedeffect.R
diff --git a/R/fitmodel.R b/R/fitmodel.R
index 24bfd0aff7c03d8694004e73eb8a6027fdff8a5d..6118363e55afa548af0faa5c47a4b4404a0a862d 100644
--- a/R/fitmodel.R
+++ b/R/fitmodel.R
@@ -13,7 +13,6 @@
 #' data(iris)
 #' formula <- Sepal.Length ~ Sepal.Width + Petal.Length
 #' isValidInput2fit(iris, formula) # Returns TRUE if all required variables are present
-#' @keywords internal
 #' @export
 isValidInput2fit <- function(data2fit, formula){
   vec_bool <- all.vars(formula) %in% colnames(data2fit)
@@ -43,7 +42,7 @@ isValidInput2fit <- function(data2fit, formula){
 #' # Drop the random effects related to 'group'
 #' modified_formula <- drop_randfx(formula)
 #'
-#' @importFrom stats terms
+#' @importFrom stats terms drop.terms
 #' @importFrom stats update
 #'
 #' @export
@@ -113,8 +112,8 @@ is_fullrank <- function(metadata, formula) {
 #' @return Fitted model object or NULL if there was an error
 #' @export
 #' @examples
-#' .fitModel(formula = mpg ~ cyl + disp, data = mtcars)
-.fitModel <- function(formula, data, ...) {
+#' fitModel(formula = mpg ~ cyl + disp, data = mtcars)
+fitModel <- function(formula, data, ...) {
   # Fit the model using glm.nb from the GLmmTMB package
   model <- glmmTMB::glmmTMB(formula, ..., data = data ) 
   model$frame <- data
@@ -140,12 +139,12 @@ is_fullrank <- function(metadata, formula) {
 #' @return Fitted model object or NULL if there was an error
 #' @export
 #' @examples
-#' .subsetData_andfit(group = "setosa", group_by = "Species", 
+#' subsetData_andfit(group = "setosa", group_by = "Species", 
 #'                  formula = Sepal.Length ~ Sepal.Width + Petal.Length, 
 #'                  data = iris )
-.subsetData_andfit <- function(group, group_by, formula, data, ...) {
+subsetData_andfit <- function(group, group_by, formula, data, ...) {
   subset_data <- data[data[[group_by]] == group, ]
-  fit_res <- .fitModel(formula, subset_data, ...)
+  fit_res <- fitModel(formula, subset_data, ...)
   #glance_df <- glance.negbin(group_by ,group , fit_res)
   #tidy_df <- tidy.negbin(group_by ,group,fit_res )
   #list(glance = glance_df, summary = tidy_df)
@@ -174,7 +173,7 @@ launchFit <- function(group, group_by, formula, data, ...) {
   tryCatch(
     expr = {
       withCallingHandlers(
-          .subsetData_andfit(group, group_by, formula, data, ...),
+          subsetData_andfit(group, group_by, formula, data, ...),
           warning = function(w) {
             message(paste(Sys.time(), "warning for group", group, ":", conditionMessage(w)))
             invokeRestart("muffleWarning")
@@ -204,14 +203,14 @@ launchFit <- function(group, group_by, formula, data, ...) {
 #' @importFrom stats setNames
 #' @export
 #' @examples
-#' .parallel_fit(group_by = "Species", "setosa", 
+#' parallel_fit(group_by = "Species", "setosa", 
 #'                formula = Sepal.Length ~ Sepal.Width + Petal.Length, 
 #'                data = iris, n.cores = 1, log_file = "log.txt" )
-.parallel_fit <- function(groups, group_by, formula, data, n.cores = NULL, log_file,  ...) {
+parallel_fit <- function(groups, group_by, formula, data, n.cores = NULL, log_file,  ...) {
   if (is.null(n.cores)) n.cores <- parallel::detectCores()
   
   clust <- parallel::makeCluster(n.cores, outfile = log_file)
-  parallel::clusterExport(clust, c(".subsetData_andfit", ".fitModel"),  envir=environment())
+  parallel::clusterExport(clust, c("subsetData_andfit", "fitModel"),  envir=environment())
   results_fit <- parallel::parLapply(clust, X = stats::setNames(groups, groups), fun = launchFit, 
                       group_by = group_by, formula = formula, data = data, ...)
                                      
@@ -244,7 +243,7 @@ fitModelParallel <- function(formula, data, group_by, n.cores = NULL, log_file =
   
   groups <- unique(data[[group_by]])
   # Fit models in parallel and capture the results
-  results <- .parallel_fit(groups, group_by, formula, data, n.cores, log_file, ...)
+  results <- parallel_fit(groups, group_by, formula, data, n.cores, log_file, ...)
   #results <- mergeListDataframes(results)
   return(results)
 }
diff --git a/R/glance_tmb.R b/R/glance_glmmtmb.R
similarity index 84%
rename from R/glance_tmb.R
rename to R/glance_glmmtmb.R
index a9bd29e16e0ad99cc675a9af7b3f6b20c5143855..1fbea83c67210c30b656343e58a20278f45417ba 100644
--- a/R/glance_tmb.R
+++ b/R/glance_glmmtmb.R
@@ -6,7 +6,7 @@
 #' This function takes a list of glmmTMB models and extracts the summary statistics (AIC, BIC, logLik, deviance,
 #' df.resid, and dispersion) for each model and returns them as a single DataFrame.
 #'
-#' @param l_tmb A list of glmmTMB models or a unique glmmTMB obj model
+#' @param list_tmb A list of glmmTMB models or a unique glmmTMB obj model
 #' @return A DataFrame with the summary statistics for all the glmmTMB models in the list.
 #' @export
 #' @importFrom stats setNames
@@ -15,10 +15,10 @@
 #' models <-  fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length, 
 #'                            group_by = "Species",n.cores = 1, data = iris)
 #' result <- glance_tmb(models)
-glance_tmb <- function(l_tmb){
-  if (identical(class(l_tmb), "glmmTMB")) return(getGlance(l_tmb))
-  l_group <- attributes(l_tmb)$names
-  l_glance <- lapply(stats::setNames(l_group, l_group), function(group) getGlance(l_tmb[[group]]))
+glance_tmb <- function(list_tmb){
+  if (identical(class(list_tmb), "glmmTMB")) return(getGlance(list_tmb))
+  l_group <- attributes(list_tmb)$names
+  l_glance <- lapply(stats::setNames(l_group, l_group), function(group) getGlance(list_tmb[[group]]))
   return(do.call("rbind", l_glance))
 }
 
diff --git a/R/identityplot.R b/R/identity_plot.R
similarity index 100%
rename from R/identityplot.R
rename to R/identity_plot.R
diff --git a/R/mock-rnaseq.R b/R/mock_rnaseq.R
similarity index 95%
rename from R/mock-rnaseq.R
rename to R/mock_rnaseq.R
index 0beccbcf91a6c2cb7b32116ca23784f9e076c4f0..6e3aef3d6f2431986810830fd7ff009fbcf72f03 100644
--- a/R/mock-rnaseq.R
+++ b/R/mock_rnaseq.R
@@ -12,8 +12,8 @@
 #' @examples
 #' matx_dispersion <- matrix(1:12, nrow = 3, ncol = 4)
 #' matx_bool_replication <- matrix(TRUE, nrow = 3, ncol = 4)
-#' .isDispersionMatrixValid(matx_dispersion, matx_bool_replication)
-.isDispersionMatrixValid <- function(matx_dispersion, matx_bool_replication) {
+#' is_dispersionMatrixValid(matx_dispersion, matx_bool_replication)
+is_dispersionMatrixValid <- function(matx_dispersion, matx_bool_replication) {
   expected_nb_column <- dim(matx_bool_replication)[2]
   if (expected_nb_column != dim(matx_dispersion)[2]) {
     return(FALSE)
@@ -81,7 +81,7 @@ mock_rnaseq <- function(list_var, n_genes, min_replicates, max_replicates, seque
   matx_Muij <- getMu_ij_matrix(df_inputSimulation)
   l_sampleID <- getSampleID(list_var)
   matx_bool_replication <- generateReplicationMatrix(list_var, min_replicates, max_replicates)
-  mu_ij_matx_rep <- .replicateMatrix(matx_Muij, matx_bool_replication)
+  mu_ij_matx_rep <- replicateMatrix(matx_Muij, matx_bool_replication)
   
   
   dispersion <- getValidDispersion(dispersion)
@@ -92,7 +92,7 @@ mock_rnaseq <- function(list_var, n_genes, min_replicates, max_replicates, seque
   
   ## same order as mu_ij_matx_rep
   matx_dispersion <- matx_dispersion[ order(row.names(matx_dispersion)), ]
-  matx_dispersion_rep <- .replicateMatrix(matx_dispersion, matx_bool_replication)
+  matx_dispersion_rep <- replicateMatrix(matx_dispersion, matx_bool_replication)
   matx_countsTable <- generateCountTable(mu_ij_matx_rep, matx_dispersion_rep)
 
   message("Counts simulation: Done")
@@ -196,8 +196,8 @@ generateReplicationMatrix <- function(list_var, min_replicates, max_replicates)
 #' @examples
 #' matrix <- matrix(1:9, nrow = 3, ncol = 3)
 #' replication_matrix <- matrix(TRUE, nrow = 3, ncol = 3)
-#' .replicateMatrix(matrix, replication_matrix)
-.replicateMatrix <- function(matrix, replication_matrix) {
+#' replicateMatrix(matrix, replication_matrix)
+replicateMatrix <- function(matrix, replication_matrix) {
   n_genes <- dim(matrix)[1]
   rep_list <- colSums(replication_matrix)
   replicated_indices <- rep(seq_len(ncol(matrix)), times = rep_list)
diff --git a/R/plot_metrics.R b/R/plot_metrics.R
index de36e10ce75da342ee5723030daf8a482742c89e..24fbb45bfc7d06bc87f5724cae3ae9b59b42b841 100644
--- a/R/plot_metrics.R
+++ b/R/plot_metrics.R
@@ -35,7 +35,7 @@ subset_glance <- function(glance_df, focus){
 #' This function generates a density plot of the specified metrics for the given
 #' list of generalized linear mixed models (GLMMs).
 #'
-#' @param l_tmb A list of GLMM objects to extract metrics from.
+#' @param list_tmb A list of GLMM objects to extract metrics from.
 #' @param focus A character vector specifying the metrics to focus on. Possible
 #'   values include "AIC", "BIC", "logLik", "deviance", "df.resid", and
 #'   "dispersion". If \code{NULL}, all available metrics will be plotted.
@@ -51,8 +51,8 @@ subset_glance <- function(glance_df, focus){
 #' models_list <-  fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length, 
 #'                      group_by = "Species",n.cores = 1, data = iris)
 #' metrics_plot(models_list, focus = c("AIC", "BIC", "deviance"))
-metrics_plot <- function(l_tmb, focus = NULL) {
-  glance_df <- glance_tmb(l_tmb)
+metrics_plot <- function(list_tmb, focus = NULL) {
+  glance_df <- glance_tmb(list_tmb)
   glance_df$group_id <- rownames(glance_df)
   if (!is.null(focus)) {
     glance_df <- subset_glance(glance_df, focus)
diff --git a/R/prepare_data2fit.R b/R/prepare_data2fit.R
index 33850af905e5cd0133408ccf91f93005434c6827..e88adf9c3575fdb2281105748ad400ef2c827fd1 100644
--- a/R/prepare_data2fit.R
+++ b/R/prepare_data2fit.R
@@ -36,8 +36,8 @@ countMatrix_2longDtf <- function(countMatrix, value_name = "kij", id_vars = "gen
 #' list_var <- init_variable()
 #' mock_data <- mock_rnaseq(list_var, n_genes = 3, 2,2, 2)
 #' dtf_countLong <- countMatrix_2longDtf(mock_data$counts)
-#' .getColumnWithSampleID(dtf_countLong, mock_data$metadata)
-.getColumnWithSampleID <- function(dtf_countsLong, metadata) {
+#' getColumnWithSampleID(dtf_countLong, mock_data$metadata)
+getColumnWithSampleID <- function(dtf_countsLong, metadata) {
   example_spleID <- as.character(dtf_countsLong[1, "sampleID"])
   regex <- paste("^", as.character(dtf_countsLong[1, "sampleID"]), "$", sep = "")
   for (indice_col in dim(metadata)[2]) {
@@ -76,7 +76,7 @@ prepareData2fit <- function(countMatrix, metadata, normalization = TRUE , respon
   }
 
   dtf_countsLong <- countMatrix_2longDtf(countMatrix, response_name)
-  metadata_columnForjoining <- .getColumnWithSampleID(dtf_countsLong, metadata)
+  metadata_columnForjoining <- getColumnWithSampleID(dtf_countsLong, metadata)
   if (is.na(metadata_columnForjoining)) {
     stop("SampleIDs do not seem to correspond between countMatrix and metadata")
   }
diff --git a/R/rocplot.R b/R/roc_plot.R
similarity index 100%
rename from R/rocplot.R
rename to R/roc_plot.R
diff --git a/R/scalingsequencingdepth.R b/R/sequencing_depth_scaling.R
similarity index 100%
rename from R/scalingsequencingdepth.R
rename to R/sequencing_depth_scaling.R
diff --git a/R/simulation.R b/R/simulation.R
index b1a736a3ae64f5fe2b4431399c8332a46d75dffc..4399ec472fbcf2e20e4411214d341efc1e433b05 100644
--- a/R/simulation.R
+++ b/R/simulation.R
@@ -73,6 +73,10 @@ getCoefficients <- function(list_var, l_dataFromMvrnorm, l_dataFromUser, n_genes
 #' @param dtf_coef The coefficient data frame.
 #' @return The coefficient data frame with log_qij column added.
 #' @export
+#' @examples
+#' list_var <- init_variable()
+#' dtf_coef <- getInput2simulation(list_var, 10)
+#' dtf_coef <- getLog_qij(dtf_coef)
 getLog_qij <- function(dtf_coef) {
   dtf_beta_numeric <- dtf_coef[sapply(dtf_coef, is.numeric)]
   dtf_coef$log_qij <- rowSums(dtf_beta_numeric, na.rm = TRUE)
@@ -114,6 +118,13 @@ getMu_ij <- function(dtf_coef) {
 
 #' @export
 #' @return A Mu_ij matrix.
+#' @examples
+#' list_var <- init_variable()
+#' dtf_coef <- getInput2simulation(list_var, 10)
+#' dtf_coef <- getLog_qij(dtf_coef)
+#' dtf_coef <- addBasalExpression(dtf_coef, 10, c(10, 20, 0))
+#' dtf_coef<- getMu_ij(dtf_coef)
+#' getMu_ij_matrix(dtf_coef)
 getMu_ij_matrix <- function(dtf_coef) {
   column_names <- colnames(dtf_coef)
   idx_var <- grepl(pattern = "label", column_names)
@@ -167,4 +178,171 @@ getSubCountsTable <- function(matx_Muij, matx_dispersion, replicateID, l_bool_re
   return(matx_kij)
 }
 
+#' getReplicationMatrix
+#'
+#' @param minN Minimum number of replicates for each sample
+#' @param maxN Maximum number of replicates for each sample
+#' @param n_samples Number of samples
+#' @export
+#' @return A replication matrix indicating which samples are replicated
+getReplicationMatrix <- function(minN, maxN, n_samples) {
+  
+  # Create a list of logical vectors representing the minimum number of replicates
+  l_replication_minimum = lapply(1:n_samples, 
+                                 FUN = function(i) rep(TRUE, times = minN) )
+  
+  # Create a list of random logical vectors representing additional replicates
+  l_replication_random = lapply(1:n_samples, 
+                                FUN = function(i) sample(x = c(TRUE, FALSE), size = maxN-minN, replace = T) )
+  
+  # Combine the replication vectors into matrices
+  matx_replication_minimum <- do.call(cbind, l_replication_minimum)
+  matx_replication_random <- do.call(cbind, l_replication_random)
+  
+  # Combine the minimum replicates and random replicates into a single matrix
+  matx_replication <- rbind(matx_replication_minimum, matx_replication_random)
+  
+  # Sort the columns of the replication matrix in descending order
+  matx_replication = apply(matx_replication, 2, sort, decreasing = TRUE ) %>% matrix(nrow = maxN)
+  
+  return(matx_replication)
+}
+
+#' getCountsTable
+#'
+#' @param matx_Muij Matrix of mean expression values for each gene and sample
+#' @param matx_dispersion Matrix of dispersion values for each gene and sample
+#' @param matx_bool_replication Replication matrix indicating which samples are replicated
+#'
+#' @return A counts table containing simulated read counts for each gene and sample
+getCountsTable <- function(matx_Muij ,  matx_dispersion, matx_bool_replication ){
+  max_replicates <-  dim(matx_bool_replication)[1]
+  
+  # Apply the getSubCountsTable function to each row of the replication matrix
+  l_countsTable = lapply(1:max_replicates, function(i) getSubCountsTable(matx_Muij , matx_dispersion, i, matx_bool_replication[i,]  ))
+  
+  # Combine the counts tables into a single matrix
+  countsTable = do.call(cbind, l_countsTable)
+  
+  return(countsTable %>% as.data.frame())
+}
+
+#' getDispersionMatrix
+#'
+#' @param list_var A list of variables (already initialized)
+#' @param n_genes Number of genes
+#' @param dispersion Vector of dispersion values for each gene
+#' @export
+#'
+#' @return A matrix of dispersion values for each gene and sample
+getDispersionMatrix <- function(list_var, n_genes, dispersion = stats::runif(n_genes, min = 0, max = 1000)){
+  l_geneID = paste("gene", 1:n_genes, sep = "")
+  l_sampleID = getSampleID(list_var) 
+  n_samples = length(l_sampleID) 
+  l_dispersion <- dispersion
+  
+  # Create a data frame for the dispersion values
+  dtf_dispersion = list(dispersion =  l_dispersion) %>% as.data.frame()
+  dtf_dispersion <- dtf_dispersion[, rep("dispersion", n_samples)]
+  rownames(dtf_dispersion) = l_geneID
+  colnames(dtf_dispersion) = l_sampleID
+  
+  matx_dispersion = dtf_dispersion %>% as.matrix()
+  
+  return(matx_dispersion)
+}
+
+
+
+
+
+#' Replicate rows of a data frame by group
+#'
+#' Replicates the rows of a data frame based on a grouping variable and replication counts for each group.
+#'
+#' @param df Data frame to replicate
+#' @param group_var Name of the grouping variable in the data frame
+#' @param rep_list Vector of replication counts for each group
+#' @return Data frame with replicated rows
+#' @examples
+#' df <- data.frame(group = c("A", "B"), value = c(1, 2))
+#' replicateByGroup(df, "group", c(2, 3))
+#'
+#' @export
+replicateByGroup <- function(df, group_var, rep_list) {
+  l_group_var <- df[[group_var]]
+  group_levels <- unique(l_group_var)
+  names(rep_list) <- group_levels
+  group_indices <- rep_list[l_group_var]
+  replicated_indices <- rep(seq_len(nrow(df)), times = group_indices)
+  replicated_df <- df[replicated_indices, ]
+  suffix_sampleID <- sequence(group_indices)
+  replicated_df[["sampleID"]] <- paste(replicated_df[["sampleID"]], suffix_sampleID, sep = "_")
+  rownames(replicated_df) <- NULL
+  return(replicated_df)
+}
+
+
+
+#' Replicate rows of a data frame
+#'
+#' Replicates the rows of a data frame by a specified factor.
+#'
+#' @param df Data frame to replicate
+#' @param n Replication factor for each row
+#' @return Data frame with replicated rows
+#' @export
+#' @examples
+#' df <- data.frame(a = 1:3, b = letters[1:3])
+#' replicateRows(df, 2)
+#'
+replicateRows <- function(df, n) {
+  indices <- rep(seq_len(nrow(df)), each = n)
+  replicated_df <- df[indices, , drop = FALSE]
+  rownames(replicated_df) <- NULL
+  return(replicated_df)
+}
+
+#' Get sample metadata
+#'
+#' Generates sample metadata based on the input variables, replication matrix, and number of genes.
+#'
+#' @param list_var A list of variables (already initialized)
+#' @param replicationMatrix Replication matrix
+#' @param n_genes Number of genes
+#' @return Data frame of sample metadata
+#' @importFrom data.table setorderv
+#' @export
+#' @examples
+#' list_var <- init_variable()
+#' n_genes <- 10
+#' replicationMatrix <- generateReplicationMatrix(list_var ,2, 3)
+#' getSampleMetadata(list_var, n_genes,  replicationMatrix)
+getSampleMetadata <- function(list_var, n_genes, replicationMatrix) {
+  l_sampleIDs = getSampleID(list_var)
+  metaData <- generateGridCombination_fromListVar(list_var)
+  metaData[] <- lapply(metaData, as.character) ## before reordering
+  data.table::setorderv(metaData, cols = colnames(metaData))
+  metaData[] <- lapply(metaData, as.factor)
+  metaData$sampleID <- l_sampleIDs
+  rep_list <- colSums(replicationMatrix)
+  metaData$sampleID <- as.character(metaData$sampleID) ## before replicating
+  sampleMetadata <- replicateByGroup(metaData, "sampleID", rep_list)
+  colnames(sampleMetadata) <- gsub("label_", "", colnames(sampleMetadata))
+  return(sampleMetadata)
+}
+
+
+#' getSampleID
+#'
+#' @param list_var A list of variables (already initialized)
+#' @export
+#' @return A sorted vector of sample IDs
+#' @examples
+#' getSampleID(init_variable())
+getSampleID <- function(list_var){
+  gridCombination <- generateGridCombination_fromListVar(list_var)
+  l_sampleID <- apply( gridCombination , 1 , paste , collapse = "_" ) %>% unname()
+  return(sort(l_sampleID))
+}
 
diff --git a/R/simulation2.R b/R/simulation2.R
deleted file mode 100644
index 7b5abfc8562e44c9c69825ea0ae5f1ef8ee7c617..0000000000000000000000000000000000000000
--- a/R/simulation2.R
+++ /dev/null
@@ -1,170 +0,0 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
-
-
-#' getReplicationMatrix
-#'
-#' @param minN Minimum number of replicates for each sample
-#' @param maxN Maximum number of replicates for each sample
-#' @param n_samples Number of samples
-#' @export
-#' @return A replication matrix indicating which samples are replicated
-getReplicationMatrix <- function(minN, maxN, n_samples) {
-  
-  # Create a list of logical vectors representing the minimum number of replicates
-  l_replication_minimum = lapply(1:n_samples, 
-                                 FUN = function(i) rep(TRUE, times = minN) )
-  
-  # Create a list of random logical vectors representing additional replicates
-  l_replication_random = lapply(1:n_samples, 
-                                FUN = function(i) sample(x = c(TRUE, FALSE), size = maxN-minN, replace = T) )
-  
-  # Combine the replication vectors into matrices
-  matx_replication_minimum <- do.call(cbind, l_replication_minimum)
-  matx_replication_random <- do.call(cbind, l_replication_random)
-  
-  # Combine the minimum replicates and random replicates into a single matrix
-  matx_replication <- rbind(matx_replication_minimum, matx_replication_random)
-  
-  # Sort the columns of the replication matrix in descending order
-  matx_replication = apply(matx_replication, 2, sort, decreasing = TRUE ) %>% matrix(nrow = maxN)
-  
-  return(matx_replication)
-}
-
-#' getCountsTable
-#'
-#' @param matx_Muij Matrix of mean expression values for each gene and sample
-#' @param matx_dispersion Matrix of dispersion values for each gene and sample
-#' @param matx_bool_replication Replication matrix indicating which samples are replicated
-#'
-#' @return A counts table containing simulated read counts for each gene and sample
-getCountsTable <- function(matx_Muij ,  matx_dispersion, matx_bool_replication ){
-  max_replicates <-  dim(matx_bool_replication)[1]
-  
-  # Apply the getSubCountsTable function to each row of the replication matrix
-  l_countsTable = lapply(1:max_replicates, function(i) getSubCountsTable(matx_Muij , matx_dispersion, i, matx_bool_replication[i,]  ))
-  
-  # Combine the counts tables into a single matrix
-  countsTable = do.call(cbind, l_countsTable)
-  
-  return(countsTable %>% as.data.frame())
-}
-
-#' getDispersionMatrix
-#'
-#' @param list_var A list of variables (already initialized)
-#' @param n_genes Number of genes
-#' @param dispersion Vector of dispersion values for each gene
-#' @export
-#'
-#' @return A matrix of dispersion values for each gene and sample
-getDispersionMatrix <- function(list_var, n_genes, dispersion = stats::runif(n_genes, min = 0, max = 1000)){
-  l_geneID = base::paste("gene", 1:n_genes, sep = "")
-  l_sampleID = getSampleID(list_var) 
-  n_samples = length(l_sampleID) 
-  l_dispersion <- dispersion
-  
-  # Create a data frame for the dispersion values
-  dtf_dispersion = list(dispersion =  l_dispersion) %>% as.data.frame()
-  dtf_dispersion <- dtf_dispersion[, rep("dispersion", n_samples)]
-  rownames(dtf_dispersion) = l_geneID
-  colnames(dtf_dispersion) = l_sampleID
-  
-  matx_dispersion = dtf_dispersion %>% as.matrix()
-  
-  return(matx_dispersion)
-}
-
-
-
-
-
-#' Replicate rows of a data frame by group
-#'
-#' Replicates the rows of a data frame based on a grouping variable and replication counts for each group.
-#'
-#' @param df Data frame to replicate
-#' @param group_var Name of the grouping variable in the data frame
-#' @param rep_list Vector of replication counts for each group
-#' @return Data frame with replicated rows
-#' @examples
-#' df <- data.frame(group = c("A", "B"), value = c(1, 2))
-#' .replicateByGroup(df, "group", c(2, 3))
-#'
-#' @export
-.replicateByGroup <- function(df, group_var, rep_list) {
-  l_group_var <- df[[group_var]]
-  group_levels <- unique(l_group_var)
-  names(rep_list) <- group_levels
-  group_indices <- rep_list[l_group_var]
-  replicated_indices <- rep(seq_len(nrow(df)), times = group_indices)
-  replicated_df <- df[replicated_indices, ]
-  suffix_sampleID <- sequence(group_indices)
-  replicated_df[["sampleID"]] <- paste(replicated_df[["sampleID"]], suffix_sampleID, sep = "_")
-  rownames(replicated_df) <- NULL
-  return(replicated_df)
-}
-
-
-
-#' Replicate rows of a data frame
-#'
-#' Replicates the rows of a data frame by a specified factor.
-#'
-#' @param df Data frame to replicate
-#' @param n Replication factor for each row
-#' @return Data frame with replicated rows
-#' @export
-#' @examples
-#' df <- data.frame(a = 1:3, b = letters[1:3])
-#' .replicateRows(df, 2)
-#'
-.replicateRows <- function(df, n) {
-  indices <- rep(seq_len(nrow(df)), each = n)
-  replicated_df <- df[indices, , drop = FALSE]
-  rownames(replicated_df) <- NULL
-  return(replicated_df)
-}
-
-#' Get sample metadata
-#'
-#' Generates sample metadata based on the input variables, replication matrix, and number of genes.
-#'
-#' @param list_var A list of variables (already initialized)
-#' @param replicationMatrix Replication matrix
-#' @param n_genes Number of genes
-#' @return Data frame of sample metadata
-#' @importFrom data.table setorderv
-#' @export
-#' @examples
-#' list_var <- init_variable()
-#' n_genes <- 10
-#' replicationMatrix <- generateReplicationMatrix(list_var ,2, 3)
-#' getSampleMetadata(list_var, n_genes,  replicationMatrix)
-getSampleMetadata <- function(list_var, n_genes, replicationMatrix) {
-  l_sampleIDs = getSampleID(list_var)
-  metaData <- generateGridCombination_fromListVar(list_var)
-  metaData[] <- lapply(metaData, as.character) ## before reordering
-  data.table::setorderv(metaData, cols = colnames(metaData))
-  metaData[] <- lapply(metaData, as.factor)
-  metaData$sampleID <- l_sampleIDs
-  rep_list <- colSums(replicationMatrix)
-  metaData$sampleID <- as.character(metaData$sampleID) ## before replicating
-  sampleMetadata <- .replicateByGroup(metaData, "sampleID", rep_list)
-  colnames(sampleMetadata) <- gsub("label_", "", colnames(sampleMetadata))
-  return(sampleMetadata)
-}
-
-
-#' getSampleID
-#'
-#' @param list_var A list of variables (already initialized)
-#' @export
-#' @return A sorted vector of sample IDs
-getSampleID <- function(list_var){
-  gridCombination <- generateGridCombination_fromListVar(list_var)
-  l_sampleID <- apply( gridCombination , 1 , paste , collapse = "_" ) %>% unname()
-  return(sort(l_sampleID))
-}
-
-
diff --git a/R/simulation_initialization.R b/R/simulation_initialization.R
index 9e0a068e5e651b840f9f96eb7ce7104ad1203b4e..65e22f671c4c3473307ad41c098b4e2a2f735313 100644
--- a/R/simulation_initialization.R
+++ b/R/simulation_initialization.R
@@ -279,47 +279,6 @@ getNumberOfCombinationsInInteraction <- function(list_var, between) {
   return(n_combinations)
 }
 
-#' getGridCombination
-#'
-#' Generates all possible combinations of labels.
-#'
-#' @param l_labels List of label vectors
-#'
-#' @return A data frame with all possible combinations of labels
-#' @export
-#'
-#' @examples
-#' l_labels <- list(
-#'   c("A", "B", "C"),
-#'   c("X", "Y")
-#' )
-#' getGridCombination(l_labels)
-getGridCombination <- function(l_labels) {
-  grid <- expand.grid(l_labels)
-  colnames(grid) <- paste("label", seq_along(l_labels), sep = "_")
-  return(grid)
-}
-
-
-
-#' Get grid combination from list_var
-#'
-#' @param list_var A list of variables (already initialized)
-#'
-#' @return
-#' The grid combination between variable in list_var
-#' @export
-generateGridCombination_fromListVar <- function (list_var){
-  l_levels <- getGivenAttribute(list_var, "level") %>% unlist()
-  vars <- names(l_levels)
-  l_levels <- l_levels[vars]
-  l_labels <- getLabels(l_variables2labelized = vars, l_nb_label = l_levels)
-  gridComb <- getGridCombination(l_labels)
-  colnames(gridComb) <- paste("label", vars, sep = "_")
-  return(gridComb)
-}
-
-
 #' Fill in interaction
 #'
 #' @param list_var A list of variables (already initialized)
@@ -351,42 +310,3 @@ fillInInteraction <- function(list_var, between, mu, sd, level) {
   return(sub_dtf)
 }
 
-#' Get the list of variable names
-#'
-#' @param input R list, e.g., output of init_variable
-#'
-#' @return
-#' A character vector with the names of variables
-getListVar <- function(input) attributes(input)$names
-
-#' Get a given attribute from a list of variables
-#'
-#' @param list_var A list of variables (already initialized)
-#' @param attribute A string specifying the attribute to retrieve in all occurrences of the list
-#'
-#' @return
-#' A list without NULL values
-getGivenAttribute <- function(list_var, attribute) {
-  l <- lapply(list_var, FUN = function(var) var[[attribute]]) 
-  l_withoutNull <- l[!vapply(l, is.null, logical(1))]
-  return(l_withoutNull)
-}
-
-
-#' Get labels for variables
-#'
-#' @param l_variables2labelized A list of variables
-#' @param l_nb_label A list of numeric values representing the number of levels per variable
-#'
-#' @return
-#' A list of labels per variable
-getLabels <- function(l_variables2labelized, l_nb_label) {
-  getVarNameLabel <- function(name, level) {
-    list_label <- paste(name, 1:level, sep = "")
-    return(list_label)
-  }
-  
-  listLabels <- Map(getVarNameLabel, l_variables2labelized, l_nb_label)
-  return(listLabels)
-}
-
diff --git a/R/simulationreport.R b/R/simulation_report.R
similarity index 96%
rename from R/simulationreport.R
rename to R/simulation_report.R
index b15cb24f8a3f1f8157e202c53c2da5843cb2f8f2..2b3598836d6c85c1c59d3b3efd4e351b901038b1 100644
--- a/R/simulationreport.R
+++ b/R/simulation_report.R
@@ -91,19 +91,19 @@ simulationReport <- function(mock_obj, list_tmb = NULL,
       TMB_comparison_df <- compareInferenceToExpected(tidyRes, mock_obj$groundTruth$effects, formula_used)
       TMB_comparison_df <- getLabelExpected(TMB_comparison_df, coeff_threshold, alt_hypothesis)
       TMB_comparison_df$from <- "HTRfit"
-      tmb_disp_inferred <- extractTMBDispersion(list_tmb)
+      tmb_disp_inferred <- extract_tmbDispersion(list_tmb)
       TMB_dispersion_df <- getDispersionComparison(tmb_disp_inferred, mock_data$groundTruth$gene_dispersion)
       TMB_dispersion_df$from <- 'HTRfit'
   }
   
   if (!is.null(dds_obj)){
-      deseq2_wrapped <- wrapper_DESeq2(dds, coeff_threshold, alt_hypothesis)
+      deseq2_wrapped <- wrap_dds(dds, coeff_threshold, alt_hypothesis)
       DESEQ_comparison_df <- inferenceToExpected_withFixedEff(deseq2_wrapped$fixEff, mock_obj$groundTruth$effects)
       DESEQ_comparison_df <- getLabelExpected(DESEQ_comparison_df, coeff_threshold, alt_hypothesis)
       DESEQ_comparison_df$from <- "DESeq2"
       DESEQ_comparison_df$component <- NA
       DESEQ_comparison_df$group <- NA
-      DESEQ_disp_inferred <- extractDESeqDispersion(deseq2_wrapped)
+      DESEQ_disp_inferred <- extract_ddsDispersion(deseq2_wrapped)
       DESEQ_dispersion_df <- getDispersionComparison(DESEQ_disp_inferred , mock_data$groundTruth$gene_dispersion)
       DESEQ_dispersion_df$from <- 'DESeq2'
   }
diff --git a/R/tidy_glmmtmb.R b/R/tidy_glmmtmb.R
index d262aa3dbfc986bf4a9edf23a801b0d66e5bf26f..7c2ebb87447e47eecc0a1a76373a8cf4028aea85 100644
--- a/R/tidy_glmmtmb.R
+++ b/R/tidy_glmmtmb.R
@@ -70,7 +70,7 @@ getTidyGlmmTMB <- function(glm_TMB, ID){
 #'
 #' This function takes a list of glmmTMB models and extracts a tidy summary of the fixed and random effects from each model. It then combines the results into a single data frame.
 #'
-#' @param l_tmb A list of glmmTMB model objects.
+#' @param list_tmb A list of glmmTMB model objects.
 #' @return A data frame containing a tidy summary of the fixed and random effects from all glmmTMB models in the list.
 #' @export
 #' @examples
@@ -78,13 +78,13 @@ getTidyGlmmTMB <- function(glm_TMB, ID){
 #' model2 <- glmmTMB::glmmTMB(Petal.Length ~ Sepal.Length + Sepal.Width + (1 | Species), data = iris)
 #' model_list <- list(Model1 = model1, Model2 = model2)
 #' tidy_summary <- tidy_tmb(model_list)
-tidy_tmb <- function(l_tmb){
-    if (identical(class(l_tmb), "glmmTMB")) return(getTidyGlmmTMB(l_tmb, "glmmTMB"))
-    attributes(l_tmb)$names
-    l_tidyRes <- lapply(attributes(l_tmb)$names,
+tidy_tmb <- function(list_tmb){
+    if (identical(class(list_tmb), "glmmTMB")) return(getTidyGlmmTMB(list_tmb, "glmmTMB"))
+    attributes(list_tmb)$names
+    l_tidyRes <- lapply(attributes(list_tmb)$names,
                  function(ID)
                    {
-                      glm_TMB <- l_tmb[[ID]]
+                      glm_TMB <- list_tmb[[ID]]
                       getTidyGlmmTMB(glm_TMB, ID)
                   }
                 )
@@ -117,22 +117,6 @@ build_missingColumn_with_na <- function(df, l_columns = c("effect", "component",
   return(df)
 }
 
-#' Remove Duplicated Words from Strings
-#'
-#' This function takes a vector of strings and removes duplicated words within each string.
-#'
-#' @param strings A character vector containing strings with potential duplicated words.
-#' @return A character vector with duplicated words removed from each string.
-#' @export
-#' @examples
-#'
-#' words <- c("hellohello", "worldworld", "programmingprogramming", "R isis great")
-#' cleaned_words <- removeDuplicatedWord(words)
-removeDuplicatedWord <- function(strings){
-  gsub("(.*)\\1+", "\\1", strings, perl = TRUE)
-}
-
-
 
 
 #' Convert Correlation Matrix to Data Frame
@@ -273,28 +257,3 @@ renameColumns <- function(df, old_names  = c("Estimate","Std..Error", "z.value",
   return(df)
 }
 
-
-
-#' Reorder the columns of a dataframe
-#'
-#' This function reorders the columns of a dataframe according to the specified column order.
-#'
-#' @param df The input dataframe.
-#' @param columnOrder A vector specifying the desired order of columns.
-#'
-#' @return A dataframe with columns reordered according to the specified column order.
-#' @export
-#' @examples
-#' # Example dataframe
-#' df <- data.frame(A = 1:3, B = 4:6, C = 7:9)
-#'
-#' # Define the desired column order
-#' columnOrder <- c("B", "C", "A")
-#'
-#' # Reorder the columns of the dataframe
-#' df <- reorderColumns(df, columnOrder)
-reorderColumns <- function(df, columnOrder) {
-  df <- df[, columnOrder, drop = FALSE]
-  return(df)
-}
-
diff --git a/R/updatefitmodel.R b/R/update_fittedmodel.R
similarity index 76%
rename from R/updatefitmodel.R
rename to R/update_fittedmodel.R
index 17370bd7bb2dd4e4fb6bc928e5030216555afd8a..e408b2a8f8a87e28881c649e6357fbc612fea603 100644
--- a/R/updatefitmodel.R
+++ b/R/update_fittedmodel.R
@@ -2,12 +2,12 @@
 
 
 
-#' Update GLMNB models in parallel.
+#' Update glmmTMB models in parallel.
 #'
-#' This function fits GLMNB models in parallel using multiple cores, allowing for faster computation.
+#' This function fits glmmTMB models in parallel using multiple cores, allowing for faster computation.
 #'
 #' @param formula Formula for the GLMNB model.
-#' @param l_tmb List of GLMNB objects.
+#' @param list_tmb List of glmmTMB objects.
 #' @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.
 #' @param ... Additional arguments to be passed to the glmmTMB::glmmTMB function.
@@ -22,24 +22,24 @@
 #' 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, l_tmb, n.cores = NULL, log_file = "log.txt", ...) {
+updateParallel <- function(formula, list_tmb, n.cores = NULL, log_file = "log.txt", ...) {
     
-    isValidInput2fit(l_tmb[[1]]$frame, formula)
+    isValidInput2fit(list_tmb[[1]]$frame, formula)
   
-    is_fullrank(l_tmb[[1]]$frame, formula)
+    is_fullrank(list_tmb[[1]]$frame, formula)
     
     # Fit models update in parallel and capture the results
-    results <- .parallel_update(formula, l_tmb, n.cores, log_file, ...)
+    results <- parallel_update(formula, list_tmb, n.cores, log_file, ...)
     return(results)
 }
 
 
-#' Internal function to fit GLMNB models in parallel.
+#' Internal function to fit glmmTMB models in parallel.
 #'
-#' This function is used internally by \code{\link{updateParallel}} to fit GLMNB models in parallel.
+#' This function is used internally by \code{\link{updateParallel}} to fit glmmTMB models in parallel.
 #'
 #' @param formula Formula for the GLMNB model.
-#' @param l_tmb List of GLMNB objects.
+#' @param list_tmb List of glmmTMB objects.
 #' @param n.cores Number of cores to use for parallel processing.
 #' @param log_file File path for the log output.
 #' @param ... Additional arguments to be passed to the glmmTMB::glmmTMB function.
@@ -52,13 +52,13 @@ updateParallel <- function(formula, l_tmb, n.cores = NULL, log_file = "log.txt",
 #' formula <- Sepal.Length ~ Sepal.Width + Petal.Length
 #' fitted_models <- fitModelParallel(formula, iris, group_by, n.cores = 1)
 #' new_formula <- Sepal.Length ~ Sepal.Width 
-#' results <- .parallel_update(new_formula, fitted_models, n.cores = 1)
-.parallel_update <- function(formula, l_tmb, n.cores = NULL, log_file = "log.txt",  ...) {
+#' results <- parallel_update(new_formula, fitted_models, n.cores = 1)
+parallel_update <- function(formula, list_tmb, n.cores = NULL, log_file = "log.txt",  ...) {
   if (is.null(n.cores)) n.cores <- parallel::detectCores()
   clust <- parallel::makeCluster(n.cores, outfile = log_file)
   #l_geneID <- attributes(l_tmb)$names
   parallel::clusterExport(clust, c("launchUpdate", "fitUpdate"),  envir=environment())
-  updated_res <- parallel::parLapply(clust, X = l_tmb, fun = launchUpdate , formula = formula, ...)
+  updated_res <- parallel::parLapply(clust, X = list_tmb, fun = launchUpdate , formula = formula, ...)
   parallel::stopCluster(clust)
   #closeAllConnections()
   return(updated_res)
@@ -69,7 +69,7 @@ updateParallel <- function(formula, l_tmb, n.cores = NULL, log_file = "log.txt",
 #'
 #' This function fits and updates a GLMNB model using the provided formula.
 #'
-#' @param glmnb_obj A GLMNB object to be updated.
+#' @param glmm_obj A glmmTMB object to be updated.
 #' @param formula Formula for the updated GLMNB model.
 #' @param ... Additional arguments to be passed to the glmmTMB::glmmTMB function.
 #' @export
@@ -83,9 +83,9 @@ updateParallel <- function(formula, l_tmb, n.cores = NULL, log_file = "log.txt",
 #' fitted_models <- fitModelParallel(formula, iris, group_by, n.cores = 1)
 #' new_formula <- Sepal.Length ~ Sepal.Width 
 #' updated_model <- fitUpdate(fitted_models[[1]], new_formula)
-fitUpdate <- function(glmnb_obj, formula , ...){
-  data = glmnb_obj$frame
-  resUpdt <- stats::update(glmnb_obj, formula, ...)
+fitUpdate <- function(glmm_obj, formula , ...){
+  data = glmm_obj$frame
+  resUpdt <- stats::update(glmm_obj, formula, ...)
   resUpdt$frame <- data
   ## family in ... => avoid error in future update
   additional_args <- list(...)
@@ -102,7 +102,7 @@ fitUpdate <- function(glmnb_obj, formula , ...){
 #'
 #' This function launches the update process for a GLMNB model, capturing and handling warnings and errors.
 #'
-#' @param glmnb_obj A GLMNB object to be updated.
+#' @param glmm_obj A glmmTMB object to be updated.
 #' @param formula Formula for the updated GLMNB model.
 #' @param ... Additional arguments to be passed to the glmmTMB::glmmTMB function.
 #' @export
@@ -116,12 +116,12 @@ fitUpdate <- function(glmnb_obj, formula , ...){
 #' fitted_models <- fitModelParallel(formula, iris, group_by, n.cores = 1)
 #' new_formula <- Sepal.Length ~ Sepal.Width 
 #' updated_model <- launchUpdate(fitted_models[[1]], new_formula)
-launchUpdate <- function(glmnb_obj, formula,  ...) {
-  group = deparse(substitute(glmnb_obj))
+launchUpdate <- function(glmm_obj, formula,  ...) {
+  group = deparse(substitute(glmm_obj))
   tryCatch(
     expr = {
       withCallingHandlers(
-        fitUpdate(glmnb_obj, formula, ...),
+        fitUpdate(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 76037f757518fe22e0760e60e26fe5b559fa9002..c25017f6d2660af7d0b879f64f791ab3acdb1a3d 100644
--- a/R/utils.R
+++ b/R/utils.R
@@ -60,6 +60,28 @@ clean_variable_name <- function(name){
     
 }
 
+#' Convert specified columns to factor
+#'
+#' @param data The input data frame
+#' @param columns The column names to be converted to factors
+#' @return The modified data frame with specified columns converted to factors
+#' @export
+#' @examples
+#' data <- data.frame( Category1 = c("A", "B", "A", "B"),
+#'                      Category2 = c("X", "Y", "X", "Z"),
+#'                      Value = 1:4,
+#'                      stringsAsFactors = FALSE )
+#' ## -- Convert columns to factors
+#' convert2Factor(data, columns = c("Category1", "Category2"))
+convert2Factor <- function(data, columns) {
+  if (is.character(columns)) {
+    columns <- match(columns, colnames(data))
+  }
+
+  if (length(columns) > 1) data[, columns] <- lapply(data[, columns], as.factor )
+  else data[, columns] <- as.factor(data[, columns])
+  return(data)
+}
 
 #' Get Setting Table
 #'
@@ -84,3 +106,163 @@ getSettingsTable <- function(n_genes, max_replicates, min_replicates, lib_size )
   return(settings_df)
 }
 
+
+#' Check if a matrix is positive definite
+#' This function checks whether a given matrix is positive definite, i.e., all of its eigenvalues are positive.
+#' @param mat The matrix to be checked.
+#' @return A logical value indicating whether the matrix is positive definite.
+#' @export
+#' @examples
+#' # Create a positive definite matrix
+#' mat1 <- matrix(c(4, 2, 2, 3), nrow = 2)
+#' is_positive_definite(mat1)
+#' # Expected output: TRUE
+#'
+#' # Create a non-positive definite matrix
+#' mat2 <- matrix(c(4, 2, 2, -3), nrow = 2)
+#' is_positive_definite(mat2)
+#' # Expected output: FALSE
+#'
+#' # Check an empty matrix
+#' mat3 <- matrix(nrow = 0, ncol = 0)
+#' is_positive_definite(mat3)
+#' # Expected output: TRUE
+#'
+is_positive_definite <- function(mat) {
+  if (nrow(mat) == 0 && ncol(mat) == 0) return(TRUE)
+  eigenvalues <- eigen(mat)$values
+  all(eigenvalues > 0)
+}
+
+
+
+
+#' Get the list of variable names
+#'
+#' @param list_var R list, e.g., output of init_variable
+#'
+#' @return
+#' A character vector with the names of variables
+#' @examples
+#' getListVar(init_variable())
+#' @export
+getListVar <- function(list_var) attributes(list_var)$names
+
+#' Get a given attribute from a list of variables
+#'
+#' @param list_var A list of variables (already initialized with init_variable)
+#' @param attribute A string specifying the attribute to retrieve in all occurrences of the list
+#' @export
+#' @return
+#' A list without NULL values
+#' @examples
+#' getGivenAttribute(init_variable(), "level")
+getGivenAttribute <- function(list_var, attribute) {
+  l <- lapply(list_var, FUN = function(var) var[[attribute]]) 
+  l_withoutNull <- l[!vapply(l, is.null, logical(1))]
+  return(l_withoutNull)
+}
+
+
+#' Get labels for variables
+#'
+#' @param l_variables2labelized A list of variables
+#' @param l_nb_label A list of numeric values representing the number of levels per variable
+#' @export
+#' @return
+#' A list of labels per variable
+#' 
+#' @examples
+#' labels <- getLabels(c("varA", "varB"), c(2, 3))
+getLabels <- function(l_variables2labelized, l_nb_label) {
+  getVarNameLabel <- function(name, level) {
+    list_label <- paste(name, 1:level, sep = "")
+    return(list_label)
+  }
+  
+  listLabels <- Map(getVarNameLabel, l_variables2labelized, l_nb_label)
+  return(listLabels)
+}
+
+
+#' getGridCombination
+#'
+#' Generates all possible combinations of labels.
+#'
+#' @param l_labels List of label vectors
+#'
+#' @return A data frame with all possible combinations of labels
+#' @export
+#'
+#' @examples
+#' l_labels <- list(
+#'   c("A", "B", "C"),
+#'   c("X", "Y")
+#' )
+#' getGridCombination(l_labels)
+getGridCombination <- function(l_labels) {
+  grid <- expand.grid(l_labels)
+  colnames(grid) <- paste("label", seq_along(l_labels), sep = "_")
+  return(grid)
+}
+
+
+
+#' Get grid combination from list_var
+#'
+#' @param list_var A list of variables (already initialized)
+#'
+#' @return
+#' The grid combination between variable in list_var
+#' @export
+#' @examples
+#' generateGridCombination_fromListVar(init_variable())
+generateGridCombination_fromListVar <- function (list_var){
+  l_levels <- getGivenAttribute(list_var, "level") %>% unlist()
+  vars <- names(l_levels)
+  l_levels <- l_levels[vars]
+  l_labels <- getLabels(l_variables2labelized = vars, l_nb_label = l_levels)
+  gridComb <- getGridCombination(l_labels)
+  colnames(gridComb) <- paste("label", vars, sep = "_")
+  return(gridComb)
+}
+
+#' Remove Duplicated Words from Strings
+#'
+#' This function takes a vector of strings and removes duplicated words within each string.
+#'
+#' @param strings A character vector containing strings with potential duplicated words.
+#' @return A character vector with duplicated words removed from each string.
+#' @export
+#' @examples
+#' words <- c("hellohello", "worldworld", "programmingprogramming", "R isis great")
+#' cleaned_words <- removeDuplicatedWord(words)
+removeDuplicatedWord <- function(strings){
+  gsub("(.*)\\1+", "\\1", strings, perl = TRUE)
+}
+
+
+#' Reorder the columns of a dataframe
+#'
+#' This function reorders the columns of a dataframe according to the specified column order.
+#'
+#' @param df The input dataframe.
+#' @param columnOrder A vector specifying the desired order of columns.
+#'
+#' @return A dataframe with columns reordered according to the specified column order.
+#' @export
+#' @examples
+#' # Example dataframe
+#' df <- data.frame(A = 1:3, B = 4:6, C = 7:9)
+#'
+#' # Define the desired column order
+#' columnOrder <- c("B", "C", "A")
+#'
+#' # Reorder the columns of the dataframe
+#' df <- reorderColumns(df, columnOrder)
+reorderColumns <- function(df, columnOrder) {
+  df <- df[, columnOrder, drop = FALSE]
+  return(df)
+}
+
+
diff --git a/R/wrapperdeseq2.R b/R/wrapper_dds.R
similarity index 93%
rename from R/wrapperdeseq2.R
rename to R/wrapper_dds.R
index f7b67a7e57fa258f9bebdc330c691a5696b59b9f..81b4f01cc2a07bd3d5069e2b29b2b2582607ee30 100644
--- a/R/wrapperdeseq2.R
+++ b/R/wrapper_dds.R
@@ -29,9 +29,9 @@
 #' dds <- DESeq2::DESeqDataSetFromMatrix(mock_data$counts , 
 #'                    mock_data$metadata, ~ genotype + environment)
 #' dds <- DESeq2::DESeq(dds, quiet = TRUE)
-#' result <- wrapper_DESeq2(dds, lfcThreshold = 1, altHypothesis = "greater")
+#' result <- wrap_dds(dds, lfcThreshold = 1, altHypothesis = "greater")
 #' @export
-wrapper_DESeq2 <- function(dds, lfcThreshold , altHypothesis, correction_method = "BH") {
+wrap_dds <- function(dds, lfcThreshold , altHypothesis, correction_method = "BH") {
   dds_full <- S4Vectors::mcols(dds) %>% as.data.frame()
   
   ## -- dispersion
@@ -40,7 +40,7 @@ wrapper_DESeq2 <- function(dds, lfcThreshold , altHypothesis, correction_method
   names(dispersion) <- rownames(dds_full)
 
   ## -- coeff
-  inference_df <- get_inference(dds_full, lfcThreshold, altHypothesis, correction_method)
+  inference_df <- get_inference_dds(dds_full, lfcThreshold, altHypothesis, correction_method)
   res <- list(dispersion = dispersion, fixEff = inference_df)
   return(res)
 }
@@ -61,13 +61,13 @@ wrapper_DESeq2 <- function(dds, lfcThreshold , altHypothesis, correction_method
 #' @examples
 #' \dontrun{
 #' # Example usage of the function
-#' inference_result <- get_inference(dds_full, lfcThreshold = 0.5, 
+#' inference_result <- get_inference_dds(dds_full, lfcThreshold = 0.5, 
 #'                                    altHypothesis = "greater", 
 #'                                    correction_method = "BH")
 #' }
 #' @importFrom stats p.adjust
 #' @export
-get_inference <- function(dds_full, lfcThreshold, altHypothesis, correction_method){
+get_inference_dds <- function(dds_full, lfcThreshold, altHypothesis, correction_method){
 
   ## -- build subdtf
   stdErr_df <- getSE_df(dds_full)
diff --git a/dev/flat_full.Rmd b/dev/flat_full.Rmd
index cb68a25a63a4f3615662f158c54442c0ca1445a5..e9425237c325cc654b55faf1283c5677a7561d24 100644
--- a/dev/flat_full.Rmd
+++ b/dev/flat_full.Rmd
@@ -85,6 +85,28 @@ clean_variable_name <- function(name){
     
 }
 
+#' Convert specified columns to factor
+#'
+#' @param data The input data frame
+#' @param columns The column names to be converted to factors
+#' @return The modified data frame with specified columns converted to factors
+#' @export
+#' @examples
+#' data <- data.frame( Category1 = c("A", "B", "A", "B"),
+#'                      Category2 = c("X", "Y", "X", "Z"),
+#'                      Value = 1:4,
+#'                      stringsAsFactors = FALSE )
+#' ## -- Convert columns to factors
+#' convert2Factor(data, columns = c("Category1", "Category2"))
+convert2Factor <- function(data, columns) {
+  if (is.character(columns)) {
+    columns <- match(columns, colnames(data))
+  }
+
+  if (length(columns) > 1) data[, columns] <- lapply(data[, columns], as.factor )
+  else data[, columns] <- as.factor(data[, columns])
+  return(data)
+}
 
 #' Get Setting Table
 #'
@@ -109,10 +131,170 @@ getSettingsTable <- function(n_genes, max_replicates, min_replicates, lib_size )
   return(settings_df)
 }
 
+
+#' Check if a matrix is positive definite
+#' This function checks whether a given matrix is positive definite, i.e., all of its eigenvalues are positive.
+#' @param mat The matrix to be checked.
+#' @return A logical value indicating whether the matrix is positive definite.
+#' @export
+#' @examples
+#' # Create a positive definite matrix
+#' mat1 <- matrix(c(4, 2, 2, 3), nrow = 2)
+#' is_positive_definite(mat1)
+#' # Expected output: TRUE
+#'
+#' # Create a non-positive definite matrix
+#' mat2 <- matrix(c(4, 2, 2, -3), nrow = 2)
+#' is_positive_definite(mat2)
+#' # Expected output: FALSE
+#'
+#' # Check an empty matrix
+#' mat3 <- matrix(nrow = 0, ncol = 0)
+#' is_positive_definite(mat3)
+#' # Expected output: TRUE
+#'
+is_positive_definite <- function(mat) {
+  if (nrow(mat) == 0 && ncol(mat) == 0) return(TRUE)
+  eigenvalues <- eigen(mat)$values
+  all(eigenvalues > 0)
+}
+
+
+
+
+#' Get the list of variable names
+#'
+#' @param list_var R list, e.g., output of init_variable
+#'
+#' @return
+#' A character vector with the names of variables
+#' @examples
+#' getListVar(init_variable())
+#' @export
+getListVar <- function(list_var) attributes(list_var)$names
+
+#' Get a given attribute from a list of variables
+#'
+#' @param list_var A list of variables (already initialized with init_variable)
+#' @param attribute A string specifying the attribute to retrieve in all occurrences of the list
+#' @export
+#' @return
+#' A list without NULL values
+#' @examples
+#' getGivenAttribute(init_variable(), "level")
+getGivenAttribute <- function(list_var, attribute) {
+  l <- lapply(list_var, FUN = function(var) var[[attribute]]) 
+  l_withoutNull <- l[!vapply(l, is.null, logical(1))]
+  return(l_withoutNull)
+}
+
+
+#' Get labels for variables
+#'
+#' @param l_variables2labelized A list of variables
+#' @param l_nb_label A list of numeric values representing the number of levels per variable
+#' @export
+#' @return
+#' A list of labels per variable
+#' 
+#' @examples
+#' labels <- getLabels(c("varA", "varB"), c(2, 3))
+getLabels <- function(l_variables2labelized, l_nb_label) {
+  getVarNameLabel <- function(name, level) {
+    list_label <- paste(name, 1:level, sep = "")
+    return(list_label)
+  }
+  
+  listLabels <- Map(getVarNameLabel, l_variables2labelized, l_nb_label)
+  return(listLabels)
+}
+
+
+#' getGridCombination
+#'
+#' Generates all possible combinations of labels.
+#'
+#' @param l_labels List of label vectors
+#'
+#' @return A data frame with all possible combinations of labels
+#' @export
+#'
+#' @examples
+#' l_labels <- list(
+#'   c("A", "B", "C"),
+#'   c("X", "Y")
+#' )
+#' getGridCombination(l_labels)
+getGridCombination <- function(l_labels) {
+  grid <- expand.grid(l_labels)
+  colnames(grid) <- paste("label", seq_along(l_labels), sep = "_")
+  return(grid)
+}
+
+
+
+#' Get grid combination from list_var
+#'
+#' @param list_var A list of variables (already initialized)
+#'
+#' @return
+#' The grid combination between variable in list_var
+#' @export
+#' @examples
+#' generateGridCombination_fromListVar(init_variable())
+generateGridCombination_fromListVar <- function (list_var){
+  l_levels <- getGivenAttribute(list_var, "level") %>% unlist()
+  vars <- names(l_levels)
+  l_levels <- l_levels[vars]
+  l_labels <- getLabels(l_variables2labelized = vars, l_nb_label = l_levels)
+  gridComb <- getGridCombination(l_labels)
+  colnames(gridComb) <- paste("label", vars, sep = "_")
+  return(gridComb)
+}
+
+#' Remove Duplicated Words from Strings
+#'
+#' This function takes a vector of strings and removes duplicated words within each string.
+#'
+#' @param strings A character vector containing strings with potential duplicated words.
+#' @return A character vector with duplicated words removed from each string.
+#' @export
+#' @examples
+#' words <- c("hellohello", "worldworld", "programmingprogramming", "R isis great")
+#' cleaned_words <- removeDuplicatedWord(words)
+removeDuplicatedWord <- function(strings){
+  gsub("(.*)\\1+", "\\1", strings, perl = TRUE)
+}
+
+
+#' Reorder the columns of a dataframe
+#'
+#' This function reorders the columns of a dataframe according to the specified column order.
+#'
+#' @param df The input dataframe.
+#' @param columnOrder A vector specifying the desired order of columns.
+#'
+#' @return A dataframe with columns reordered according to the specified column order.
+#' @export
+#' @examples
+#' # Example dataframe
+#' df <- data.frame(A = 1:3, B = 4:6, C = 7:9)
+#'
+#' # Define the desired column order
+#' columnOrder <- c("B", "C", "A")
+#'
+#' # Reorder the columns of the dataframe
+#' df <- reorderColumns(df, columnOrder)
+reorderColumns <- function(df, columnOrder) {
+  df <- df[, columnOrder, drop = FALSE]
+  return(df)
+}
+
+
 ```
 
 
-```{r test-dataFromUser}
+```{r test-utils}
 # Test unitaires pour la fonction join_dtf
 test_that("join_dtf réalise la jointure correctement", {
   # Création de données de test
@@ -142,6 +324,90 @@ test_that("clean_variable_name handles reserved names properly", {
   expect_error(clean_variable_name("interactions"))
   expect_error(clean_variable_name("correlations"))
 })
+
+
+test_that("getLabels generates labels for variables", {
+  labels <- getLabels(c("varA", "varB"), c(2, 3))
+  expect_equal(length(labels), 2)
+  expect_equal(length(labels[[1]]), 2)
+  expect_equal(length(labels[[2]]), 3)
+})
+
+test_that("getGridCombination generates a grid of combinations", {
+  labels <- list(A = c("A1", "A2"), B = c("B1", "B2", "B3"))
+  grid_combination <- getGridCombination(labels)
+  expect_equal(dim(grid_combination), c(6, 2))
+})
+
+
+test_that("generateGridCombination_fromListVar returns expected output", {
+  result <- generateGridCombination_fromListVar(init_variable())
+  expect <- data.frame(label_myVariable = c("myVariable1", "myVariable2"))
+  expect_equal(nrow(result), nrow(expect))
+  expect_equal(ncol(result), ncol(expect))
+})
+
+# Tests for convert2Factor
+test_that("convert2Factor converts specified columns to factors", {
+  # Create a sample data frame
+  data <- data.frame(
+    Category1 = c("A", "B", "A", "B"),
+    Category2 = c("X", "Y", "X", "Z"),
+    Value = 1:4,
+    stringsAsFactors = FALSE
+  )
+  
+  # Convert columns to factors
+  result <- convert2Factor(data, columns = c("Category1", "Category2"))
+  
+  # Check the output
+  expect_is(result$Category1, "factor")  # Category1 column converted to factor
+  expect_is(result$Category2, "factor")  # Category2 column converted to factor
+})
+
+test_that("removeDuplicatedWord returns expected output", {
+  words <- c("hellohello", "worldworld", "programmingprogramming", "R isis great")
+  cleaned_words <- removeDuplicatedWord(words)
+  expect_equal(cleaned_words, c("hello", "world", "programming", "R is great"))
+})
+
+test_that("reorderColumns returns expected output",{
+    df <- data.frame(A = 1:3, B = 4:6, C = 7:9)
+    # Define the desired column order
+    columnOrder <- c("B", "C", "A")
+    # Reorder the columns of the dataframe
+    df_reorder <- reorderColumns(df, columnOrder)
+    expect_equal(colnames(df_reorder), columnOrder)
+    expect_equal(dim(df_reorder), dim(df))
+
+})
+
+
+test_that( "generateGridCombination_fromListVar return expected output", {
+    ## case 1
+    gridcom <- generateGridCombination_fromListVar(init_variable())
+    expect_s3_class(gridcom, "data.frame")
+    expect_equal(gridcom$label_myVariable, factor(c("myVariable1", "myVariable2")))
+
+    ## case 2
+    init_variables <- init_variable() %>% init_variable(name = "var" , mu = 2, sd = 1, level = 3) 
+    gridcom <- generateGridCombination_fromListVar(init_variables)
+    expect_s3_class(gridcom, "data.frame")
+    expect_equal(unique(gridcom$label_myVariable), factor(c("myVariable1", "myVariable2")))
+    expect_equal(unique(gridcom$label_var), factor(c("var1", "var2", "var3")))
+
+  })
+
+test_that( "getGivenAttribute return expected output", {
+  ## -- case 1
+  level_attr <- getGivenAttribute(init_variable(), "level")
+  expect_equal(level_attr$myVariable, 2)
+
+  ## -- case 2
+  init_variables <- init_variable() %>% init_variable(name = "var" , mu = 2, sd = 1, level = 3) 
+  mu_attr <- getGivenAttribute(init_variables, "mu")
+  expect_equal(mu_attr$var, 2)
+} )
 ```
 
 
@@ -425,47 +691,6 @@ getNumberOfCombinationsInInteraction <- function(list_var, between) {
   return(n_combinations)
 }
 
-#' getGridCombination
-#'
-#' Generates all possible combinations of labels.
-#'
-#' @param l_labels List of label vectors
-#'
-#' @return A data frame with all possible combinations of labels
-#' @export
-#'
-#' @examples
-#' l_labels <- list(
-#'   c("A", "B", "C"),
-#'   c("X", "Y")
-#' )
-#' getGridCombination(l_labels)
-getGridCombination <- function(l_labels) {
-  grid <- expand.grid(l_labels)
-  colnames(grid) <- paste("label", seq_along(l_labels), sep = "_")
-  return(grid)
-}
-
-
-
-#' Get grid combination from list_var
-#'
-#' @param list_var A list of variables (already initialized)
-#'
-#' @return
-#' The grid combination between variable in list_var
-#' @export
-generateGridCombination_fromListVar <- function (list_var){
-  l_levels <- getGivenAttribute(list_var, "level") %>% unlist()
-  vars <- names(l_levels)
-  l_levels <- l_levels[vars]
-  l_labels <- getLabels(l_variables2labelized = vars, l_nb_label = l_levels)
-  gridComb <- getGridCombination(l_labels)
-  colnames(gridComb) <- paste("label", vars, sep = "_")
-  return(gridComb)
-}
-
-
 #' Fill in interaction
 #'
 #' @param list_var A list of variables (already initialized)
@@ -497,45 +722,6 @@ fillInInteraction <- function(list_var, between, mu, sd, level) {
   return(sub_dtf)
 }
 
-#' Get the list of variable names
-#'
-#' @param input R list, e.g., output of init_variable
-#'
-#' @return
-#' A character vector with the names of variables
-getListVar <- function(input) attributes(input)$names
-
-#' Get a given attribute from a list of variables
-#'
-#' @param list_var A list of variables (already initialized)
-#' @param attribute A string specifying the attribute to retrieve in all occurrences of the list
-#'
-#' @return
-#' A list without NULL values
-getGivenAttribute <- function(list_var, attribute) {
-  l <- lapply(list_var, FUN = function(var) var[[attribute]]) 
-  l_withoutNull <- l[!vapply(l, is.null, logical(1))]
-  return(l_withoutNull)
-}
-
-
-#' Get labels for variables
-#'
-#' @param l_variables2labelized A list of variables
-#' @param l_nb_label A list of numeric values representing the number of levels per variable
-#'
-#' @return
-#' A list of labels per variable
-getLabels <- function(l_variables2labelized, l_nb_label) {
-  getVarNameLabel <- function(name, level) {
-    list_label <- paste(name, 1:level, sep = "")
-    return(list_label)
-  }
-  
-  listLabels <- Map(getVarNameLabel, l_variables2labelized, l_nb_label)
-  return(listLabels)
-}
-
 ```
 
 
@@ -630,13 +816,6 @@ test_that("build_sub_obj_return_to_user returns the expected output", {
   
 })
 
-test_that("generateGridCombination_fromListVar returns expected output", {
-  result <- generateGridCombination_fromListVar(init_variable())
-  expect <- data.frame(label_myVariable = c("myVariable1", "myVariable2"))
-  expect_equal(nrow(result), nrow(expect))
-  expect_equal(ncol(result), ncol(expect))
-})
-
 test_that("add_interaction adds an interaction between variables", {
   list_var <- init_variable(name = "varA", mu = 1, sd = 1, level = 2)
   list_var <- init_variable(list_var, name = "varB", mu = 2, sd = 1, level = 3)
@@ -655,20 +834,6 @@ test_that("getNumberOfCombinationsInInteraction calculates the number of combina
   list_var <- init_variable(list_var, name = "varB", mu = 2, sd = 1, level = 3)
   expect_equal(getNumberOfCombinationsInInteraction(list_var, c("varA", "varB")), 6)
 })
-
-test_that("getLabels generates labels for variables", {
-  labels <- getLabels(c("varA", "varB"), c(2, 3))
-  expect_equal(length(labels), 2)
-  expect_equal(length(labels[[1]]), 2)
-  expect_equal(length(labels[[2]]), 3)
-})
-
-test_that("getGridCombination generates a grid of combinations", {
-  labels <- list(A = c("A1", "A2"), B = c("B1", "B2", "B3"))
-  grid_combination <- getGridCombination(labels)
-  expect_equal(dim(grid_combination), c(6, 2))
-})
-
 ```
 
 ```{r function-mvrnorm, filename = "datafrommvrnorm_manipulations" }
@@ -766,37 +931,6 @@ fillInCovarMatrice <- function(covarMatrice, covar){
   return(covarMatrice)
 }
 
-
-#' Check if a matrix is positive definite
-#' This function checks whether a given matrix is positive definite, i.e., all of its eigenvalues are positive.
-#' @param mat The matrix to be checked.
-#' @return A logical value indicating whether the matrix is positive definite.
-#' @export
-#' @examples
-#' # Create a positive definite matrix
-#' mat1 <- matrix(c(4, 2, 2, 3), nrow = 2)
-#' is_positive_definite(mat1)
-#' # Expected output: TRUE
-#'
-#' # Create a non-positive definite matrix
-#' mat2 <- matrix(c(4, 2, 2, -3), nrow = 2)
-#' is_positive_definite(mat2)
-#' # Expected output: FALSE
-#'
-#' # Check an empty matrix
-#' mat3 <- matrix(nrow = 0, ncol = 0)
-#' is_positive_definite(mat3)
-#' # Expected output: TRUE
-#'
-#' @export
-is_positive_definite <- function(mat) {
-  if (nrow(mat) == 0 && ncol(mat) == 0) return(TRUE)
-  eigenvalues <- eigen(mat)$values
-  all(eigenvalues > 0)
-}
-
-
-
 #' getGeneMetadata
 #'
 #' @inheritParams init_variable
@@ -813,7 +947,7 @@ is_positive_definite <- function(mat) {
 getGeneMetadata <- function(list_var, n_genes) {
   metaData <- generateGridCombination_fromListVar(list_var)
   n_combinations <- dim(metaData)[1]
-  genes_vec <- base::paste("gene", 1:n_genes, sep = "")
+  genes_vec <- paste("gene", 1:n_genes, sep = "")
   geneID <- rep(genes_vec, each = n_combinations)
   metaData <- cbind(geneID, metaData)
   
@@ -862,7 +996,7 @@ getDataFromMvrnorm <- function(list_var, input2mvrnorm, n_genes = 1) {
 #' samples generated from multivariate normal distribution
 #' 
 #' @export
-#'
+#' @importFrom MASS mvrnorm
 #' @examples
 #' n <- 100
 #' mu <- c(0, 0)
@@ -870,7 +1004,6 @@ getDataFromMvrnorm <- function(list_var, input2mvrnorm, n_genes = 1) {
 #' samples <- samplingFromMvrnorm(n_samplings = n, l_mu = mu, matx_cov = covMatrix)
 samplingFromMvrnorm <- function(n_samplings, l_mu, matx_cov) {
   mvrnormSamp <-  MASS::mvrnorm(n = n_samplings, mu = l_mu, Sigma = matx_cov, empirical = TRUE)
-  
   return(mvrnormSamp)
 }
 
@@ -1023,7 +1156,7 @@ test_that("getDataFromUser renvoie les données appropriées", {
 })
 ```
 
-```{r function-setCorrelation, filename =  "setCorrelation"}
+```{r function-setcorrelation, filename =  "setcorrelation"}
 
 #' Compute Covariation from Correlation and Standard Deviations
 #'
@@ -1274,6 +1407,10 @@ getCoefficients <- function(list_var, l_dataFromMvrnorm, l_dataFromUser, n_genes
 #' @param dtf_coef The coefficient data frame.
 #' @return The coefficient data frame with log_qij column added.
 #' @export
+#' @examples
+#' list_var <- init_variable()
+#' dtf_coef <- getInput2simulation(list_var, 10)
+#' dtf_coef <- getLog_qij(dtf_coef)
 getLog_qij <- function(dtf_coef) {
   dtf_beta_numeric <- dtf_coef[sapply(dtf_coef, is.numeric)]
   dtf_coef$log_qij <- rowSums(dtf_beta_numeric, na.rm = TRUE)
@@ -1315,6 +1452,13 @@ getMu_ij <- function(dtf_coef) {
 
 #' @export
 #' @return A Mu_ij matrix.
+#' @examples
+#' list_var <- init_variable()
+#' dtf_coef <- getInput2simulation(list_var, 10)
+#' dtf_coef <- getLog_qij(dtf_coef)
+#' dtf_coef <- addBasalExpression(dtf_coef, 10, c(10, 20, 0))
+#' dtf_coef<- getMu_ij(dtf_coef)
+#' getMu_ij_matrix(dtf_coef)
 getMu_ij_matrix <- function(dtf_coef) {
   column_names <- colnames(dtf_coef)
   idx_var <- grepl(pattern = "label", column_names)
@@ -1354,86 +1498,19 @@ getSubCountsTable <- function(matx_Muij, matx_dispersion, replicateID, l_bool_re
   
   if (!any(l_bool_replication))
     return(NULL) 
-  
-  matx_Muij <- matx_Muij[, l_bool_replication, drop = FALSE]
-  matx_dispersion <- matx_dispersion[, l_bool_replication, drop = FALSE] 
-  l_sampleID <- colnames(matx_Muij)
-  l_geneID <- rownames(matx_Muij)
-  dimension_mtx <- dim(matx_Muij)
-  n_genes <- dimension_mtx[1]
-  n_samples <- dimension_mtx[2]
-  matx_kij <- getKijMatrix(matx_Muij, matx_dispersion, n_genes, n_samples)
-  colnames(matx_kij) <- paste(l_sampleID, replicateID, sep = "_")
-  rownames(matx_kij) <- l_geneID
-  return(matx_kij)
-}
-
-
-```
-
-```{r test-simulation}
-
-
-# Test case 1: Check if the function returns a data frame
-test_that("getInput2simulation returns a data frame", {
-  list_var <- init_variable()
-  result <- getInput2simulation(list_var)
-  expect_is(result, "data.frame")
-  expected <- data.frame(geneID = c("gene1", "gene1"), label_myVariable = as.factor(c("myVariable1", "myVariable2")), myVariable = c(2,3))
-  expect_equal(result, expected)
-  })
-
-# Test for getCoefficients function
-test_that("getCoefficients returns the correct output", {
-  # Create dummy data
-  n_genes <- 3
-  list_var = init_variable()
-  # Call the function
-  coefficients <- getCoefficients(list_var, list(), list(), n_genes)
-  
-  # Check the output
-  expect_equal(nrow(coefficients), n_genes*list_var$myVariable$level)
-  expect_equal(colnames(coefficients), c("geneID", "label_myVariable")) 
-})
-
-# Test for getMu_ij_matrix function
-test_that("getMu_ij_matrix returns the correct output", {
-  # Create a dummy coefficients dataframe
-  dtf_coef <- data.frame(geneID = c("Gene1", "Gene1", "Gene1"),
-                         label_varA = c("A1", "A2", "A3"),
-                         label_varB = c("B1", "B2", "B3"),
-                         mu_ij = c(1, 2, 3))
-  
-  # Call the function
-  mu_matrix <- getMu_ij_matrix(dtf_coef)
-  # Check the output
-  expect_equal(dim(mu_matrix), c(1, 9)) 
-  
-})
-
-# Test for getSubCountsTable function
-test_that("getSubCountsTable returns the correct output", {
-  # Create dummy data
-  l_genes <- c("gene1", "gene2", "gene3")
-  matx_Muij = data.frame(sple1 = c(1,3,4), sple2 = c(2, 0, 9), sple3 = c(1, 69, 2)) %>% as.matrix()
-  rownames(matx_Muij) <- l_genes
-  matx_dispersion <- matrix(0.5, nrow = 3, ncol = 3)
-  replicateID <- 1
-  l_bool_replication <- c(TRUE, FALSE, TRUE)
-  
-  # Call the function
-  subcounts_table <- getSubCountsTable(matx_Muij, matx_dispersion, 1, l_bool_replication)
-  
-  # Check the output
-  expect_equal(dim(subcounts_table), c(3, 2))
-  expect_equal(rownames(subcounts_table), l_genes)
-})
-
-
-```
-
-
-```{r function-simulation2 , filename = "simulation2"}
+  
+  matx_Muij <- matx_Muij[, l_bool_replication, drop = FALSE]
+  matx_dispersion <- matx_dispersion[, l_bool_replication, drop = FALSE] 
+  l_sampleID <- colnames(matx_Muij)
+  l_geneID <- rownames(matx_Muij)
+  dimension_mtx <- dim(matx_Muij)
+  n_genes <- dimension_mtx[1]
+  n_samples <- dimension_mtx[2]
+  matx_kij <- getKijMatrix(matx_Muij, matx_dispersion, n_genes, n_samples)
+  colnames(matx_kij) <- paste(l_sampleID, replicateID, sep = "_")
+  rownames(matx_kij) <- l_geneID
+  return(matx_kij)
+}
 
 #' getReplicationMatrix
 #'
@@ -1493,7 +1570,7 @@ getCountsTable <- function(matx_Muij ,  matx_dispersion, matx_bool_replication )
 #'
 #' @return A matrix of dispersion values for each gene and sample
 getDispersionMatrix <- function(list_var, n_genes, dispersion = stats::runif(n_genes, min = 0, max = 1000)){
-  l_geneID = base::paste("gene", 1:n_genes, sep = "")
+  l_geneID = paste("gene", 1:n_genes, sep = "")
   l_sampleID = getSampleID(list_var) 
   n_samples = length(l_sampleID) 
   l_dispersion <- dispersion
@@ -1523,10 +1600,10 @@ getDispersionMatrix <- function(list_var, n_genes, dispersion = stats::runif(n_g
 #' @return Data frame with replicated rows
 #' @examples
 #' df <- data.frame(group = c("A", "B"), value = c(1, 2))
-#' .replicateByGroup(df, "group", c(2, 3))
+#' replicateByGroup(df, "group", c(2, 3))
 #'
 #' @export
-.replicateByGroup <- function(df, group_var, rep_list) {
+replicateByGroup <- function(df, group_var, rep_list) {
   l_group_var <- df[[group_var]]
   group_levels <- unique(l_group_var)
   names(rep_list) <- group_levels
@@ -1551,9 +1628,9 @@ getDispersionMatrix <- function(list_var, n_genes, dispersion = stats::runif(n_g
 #' @export
 #' @examples
 #' df <- data.frame(a = 1:3, b = letters[1:3])
-#' .replicateRows(df, 2)
+#' replicateRows(df, 2)
 #'
-.replicateRows <- function(df, n) {
+replicateRows <- function(df, n) {
   indices <- rep(seq_len(nrow(df)), each = n)
   replicated_df <- df[indices, , drop = FALSE]
   rownames(replicated_df) <- NULL
@@ -1584,7 +1661,7 @@ getSampleMetadata <- function(list_var, n_genes, replicationMatrix) {
   metaData$sampleID <- l_sampleIDs
   rep_list <- colSums(replicationMatrix)
   metaData$sampleID <- as.character(metaData$sampleID) ## before replicating
-  sampleMetadata <- .replicateByGroup(metaData, "sampleID", rep_list)
+  sampleMetadata <- replicateByGroup(metaData, "sampleID", rep_list)
   colnames(sampleMetadata) <- gsub("label_", "", colnames(sampleMetadata))
   return(sampleMetadata)
 }
@@ -1595,16 +1672,72 @@ getSampleMetadata <- function(list_var, n_genes, replicationMatrix) {
 #' @param list_var A list of variables (already initialized)
 #' @export
 #' @return A sorted vector of sample IDs
+#' @examples
+#' getSampleID(init_variable())
 getSampleID <- function(list_var){
   gridCombination <- generateGridCombination_fromListVar(list_var)
   l_sampleID <- apply( gridCombination , 1 , paste , collapse = "_" ) %>% unname()
   return(sort(l_sampleID))
 }
 
-
 ```
 
-```{r test-simulations}
+```{r test-simulation}
+# Test case 1: Check if the function returns a data frame
+test_that("getInput2simulation returns a data frame", {
+  list_var <- init_variable()
+  result <- getInput2simulation(list_var)
+  expect_is(result, "data.frame")
+  expected <- data.frame(geneID = c("gene1", "gene1"), label_myVariable = as.factor(c("myVariable1", "myVariable2")), myVariable = c(2,3))
+  expect_equal(result, expected)
+  })
+
+# Test for getCoefficients function
+test_that("getCoefficients returns the correct output", {
+  # Create dummy data
+  n_genes <- 3
+  list_var = init_variable()
+  # Call the function
+  coefficients <- getCoefficients(list_var, list(), list(), n_genes)
+  
+  # Check the output
+  expect_equal(nrow(coefficients), n_genes*list_var$myVariable$level)
+  expect_equal(colnames(coefficients), c("geneID", "label_myVariable")) 
+})
+
+# Test for getMu_ij_matrix function
+test_that("getMu_ij_matrix returns the correct output", {
+  # Create a dummy coefficients dataframe
+  dtf_coef <- data.frame(geneID = c("Gene1", "Gene1", "Gene1"),
+                         label_varA = c("A1", "A2", "A3"),
+                         label_varB = c("B1", "B2", "B3"),
+                         mu_ij = c(1, 2, 3))
+  
+  # Call the function
+  mu_matrix <- getMu_ij_matrix(dtf_coef)
+  # Check the output
+  expect_equal(dim(mu_matrix), c(1, 9)) 
+  
+})
+
+# Test for getSubCountsTable function
+test_that("getSubCountsTable returns the correct output", {
+  # Create dummy data
+  l_genes <- c("gene1", "gene2", "gene3")
+  matx_Muij = data.frame(sple1 = c(1,3,4), sple2 = c(2, 0, 9), sple3 = c(1, 69, 2)) %>% as.matrix()
+  rownames(matx_Muij) <- l_genes
+  matx_dispersion <- matrix(0.5, nrow = 3, ncol = 3)
+  replicateID <- 1
+  l_bool_replication <- c(TRUE, FALSE, TRUE)
+  
+  # Call the function
+  subcounts_table <- getSubCountsTable(matx_Muij, matx_dispersion, 1, l_bool_replication)
+  
+  # Check the output
+  expect_equal(dim(subcounts_table), c(3, 2))
+  expect_equal(rownames(subcounts_table), l_genes)
+})
+
 
 test_that("getReplicationMatrix returns the correct replication matrix", {
   minN <- 2
@@ -1709,9 +1842,9 @@ test_that("getSampleMetadata returns expected output", {
 })
 
 
-test_that(".replicateByGroup return the correct ouptut", {
+test_that("replicateByGroup return the correct ouptut", {
   df <- data.frame(group = c("A", "B"), value = c(1, 2))
-  result <- .replicateByGroup(df, "group", c(2, 3))
+  result <- replicateByGroup(df, "group", c(2, 3))
   
   expect <- data.frame(group = c("A", "A", "B", "B", "B"), 
                        value = c(1, 1, 2,2,2), 
@@ -1733,11 +1866,52 @@ test_that("getDispersionMatrix returns the correct dispersion matrix", {
 })
 
 
+# Test case: Valid input vector with numeric and positive values
+test_that("Valid input vector with numeric and positive values", {
+  input_vector <- c(0.5, 1.2, 0.8)
+  result <- getValidDispersion(input_vector)
+  expect_identical(result, input_vector)
+})
+
+# Test case: Valid input vector with numeric, positive, and negative values
+test_that("Valid input vector with numeric, positive, and negative values", {
+  input_vector <- c(0.5, -0.3, 1.2, 0.8)
+  result <- getValidDispersion(input_vector)
+  expect_identical(result, c(0.5, 1.2, 0.8))
+})
+
+# Test case: Valid input vector with non-numeric elements
+test_that("Valid input vector with non-numeric elements", {
+  input_vector <- c(0.5, "invalid", 0.8)
+  result <- getValidDispersion(input_vector)
+  expect_identical(result, c(0.5, 0.8))
+})
+
+# Test case: Empty input vector
+test_that("Empty input vector", {
+  input_vector <- numeric(0)
+  expect_error(getValidDispersion(input_vector), "Invalid dispersion values provided.")
+})
+
+# Test case: unique value in vector
+test_that("unique value in vector", {
+  input_vector <- 5
+  expect_equal(getValidDispersion(input_vector), 5)
+})
+
+# Test case: All negative values
+test_that("All negative values", {
+  input_vector <- c(-0.5, -1.2, -0.8)
+  expect_error(getValidDispersion(input_vector), "Invalid dispersion values provided.")
+})
+
+
+
 
 ```
 
 
-```{r function-mock , filename = "mock-rnaSeq" }
+```{r function-mock , filename = "mock_rnaseq" }
 
 #' Check the validity of the dispersion matrix
 #'
@@ -1750,8 +1924,8 @@ test_that("getDispersionMatrix returns the correct dispersion matrix", {
 #' @examples
 #' matx_dispersion <- matrix(1:12, nrow = 3, ncol = 4)
 #' matx_bool_replication <- matrix(TRUE, nrow = 3, ncol = 4)
-#' .isDispersionMatrixValid(matx_dispersion, matx_bool_replication)
-.isDispersionMatrixValid <- function(matx_dispersion, matx_bool_replication) {
+#' is_dispersionMatrixValid(matx_dispersion, matx_bool_replication)
+is_dispersionMatrixValid <- function(matx_dispersion, matx_bool_replication) {
   expected_nb_column <- dim(matx_bool_replication)[2]
   if (expected_nb_column != dim(matx_dispersion)[2]) {
     return(FALSE)
@@ -1819,7 +1993,7 @@ mock_rnaseq <- function(list_var, n_genes, min_replicates, max_replicates, seque
   matx_Muij <- getMu_ij_matrix(df_inputSimulation)
   l_sampleID <- getSampleID(list_var)
   matx_bool_replication <- generateReplicationMatrix(list_var, min_replicates, max_replicates)
-  mu_ij_matx_rep <- .replicateMatrix(matx_Muij, matx_bool_replication)
+  mu_ij_matx_rep <- replicateMatrix(matx_Muij, matx_bool_replication)
   
   
   dispersion <- getValidDispersion(dispersion)
@@ -1830,7 +2004,7 @@ mock_rnaseq <- function(list_var, n_genes, min_replicates, max_replicates, seque
   
   ## same order as mu_ij_matx_rep
   matx_dispersion <- matx_dispersion[ order(row.names(matx_dispersion)), ]
-  matx_dispersion_rep <- .replicateMatrix(matx_dispersion, matx_bool_replication)
+  matx_dispersion_rep <- replicateMatrix(matx_dispersion, matx_bool_replication)
   matx_countsTable <- generateCountTable(mu_ij_matx_rep, matx_dispersion_rep)
 
   message("Counts simulation: Done")
@@ -1934,8 +2108,8 @@ generateReplicationMatrix <- function(list_var, min_replicates, max_replicates)
 #' @examples
 #' matrix <- matrix(1:9, nrow = 3, ncol = 3)
 #' replication_matrix <- matrix(TRUE, nrow = 3, ncol = 3)
-#' .replicateMatrix(matrix, replication_matrix)
-.replicateMatrix <- function(matrix, replication_matrix) {
+#' replicateMatrix(matrix, replication_matrix)
+replicateMatrix <- function(matrix, replication_matrix) {
   n_genes <- dim(matrix)[1]
   rep_list <- colSums(replication_matrix)
   replicated_indices <- rep(seq_len(ncol(matrix)), times = rep_list)
@@ -1948,59 +2122,19 @@ generateReplicationMatrix <- function(list_var, min_replicates, max_replicates)
 
 ```
 
-```{r test-hiddenFunction}
-
-# Test case: Valid input vector with numeric and positive values
-test_that("Valid input vector with numeric and positive values", {
-  input_vector <- c(0.5, 1.2, 0.8)
-  result <- getValidDispersion(input_vector)
-  expect_identical(result, input_vector)
-})
-
-# Test case: Valid input vector with numeric, positive, and negative values
-test_that("Valid input vector with numeric, positive, and negative values", {
-  input_vector <- c(0.5, -0.3, 1.2, 0.8)
-  result <- getValidDispersion(input_vector)
-  expect_identical(result, c(0.5, 1.2, 0.8))
-})
-
-# Test case: Valid input vector with non-numeric elements
-test_that("Valid input vector with non-numeric elements", {
-  input_vector <- c(0.5, "invalid", 0.8)
-  result <- getValidDispersion(input_vector)
-  expect_identical(result, c(0.5, 0.8))
-})
-
-# Test case: Empty input vector
-test_that("Empty input vector", {
-  input_vector <- numeric(0)
-  expect_error(getValidDispersion(input_vector), "Invalid dispersion values provided.")
-})
-
-# Test case: unique value in vector
-test_that("unique value in vector", {
-  input_vector <- 5
-  expect_equal(getValidDispersion(input_vector), 5)
-})
-
-# Test case: All negative values
-test_that("All negative values", {
-  input_vector <- c(-0.5, -1.2, -0.8)
-  expect_error(getValidDispersion(input_vector), "Invalid dispersion values provided.")
-})
+```{r test-mock}
 
-
-# Test for .isDispersionMatrixValid
-test_that(".isDispersionMatrixValid returns TRUE for valid dimensions", {
+# Test for is_dispersionMatrixValid
+test_that("is_dispersionMatrixValid returns TRUE for valid dimensions", {
   matx_dispersion <- matrix(1:6, nrow = 2, ncol = 3)
   matx_bool_replication <- matrix(TRUE, nrow = 2, ncol = 3)
-  expect_true(.isDispersionMatrixValid(matx_dispersion, matx_bool_replication))
+  expect_true(is_dispersionMatrixValid(matx_dispersion, matx_bool_replication))
 })
 
-test_that(".isDispersionMatrixValid throws an error for invalid dimensions", {
+test_that("is_dispersionMatrixValid throws an error for invalid dimensions", {
   matx_dispersion <- matrix(1:4, nrow = 2, ncol = 2)
   matx_bool_replication <- matrix(TRUE, nrow = 2, ncol = 3)
-  expect_false(.isDispersionMatrixValid(matx_dispersion, matx_bool_replication))
+  expect_false(is_dispersionMatrixValid(matx_dispersion, matx_bool_replication))
 })
 
 # Test for generateCountTable
@@ -2013,17 +2147,14 @@ test_that("generateCountTable generates count table with correct dimensions", {
 
 
 
-# Test for .replicateMatrix
-test_that(".replicateMatrix replicates matrix correctly", {
+# Test for replicateMatrix
+test_that("replicateMatrix replicates matrix correctly", {
   matrix <- matrix(1:9, nrow = 3, ncol = 3)
   replication_matrix <- matrix(TRUE, nrow = 3, ncol = 3)
-  replicated_matrix <- .replicateMatrix(matrix, replication_matrix)
+  replicated_matrix <- replicateMatrix(matrix, replication_matrix)
   expect_equal(dim(replicated_matrix), c(3, 9))
 })
 
-```
-
-```{r  test-mock}
 
 # Test for mock_rnaseq
 #test_that("mock_rnaseq returns expected output", {
@@ -2100,8 +2231,8 @@ countMatrix_2longDtf <- function(countMatrix, value_name = "kij", id_vars = "gen
 #' list_var <- init_variable()
 #' mock_data <- mock_rnaseq(list_var, n_genes = 3, 2,2, 2)
 #' dtf_countLong <- countMatrix_2longDtf(mock_data$counts)
-#' .getColumnWithSampleID(dtf_countLong, mock_data$metadata)
-.getColumnWithSampleID <- function(dtf_countsLong, metadata) {
+#' getColumnWithSampleID(dtf_countLong, mock_data$metadata)
+getColumnWithSampleID <- function(dtf_countsLong, metadata) {
   example_spleID <- as.character(dtf_countsLong[1, "sampleID"])
   regex <- paste("^", as.character(dtf_countsLong[1, "sampleID"]), "$", sep = "")
   for (indice_col in dim(metadata)[2]) {
@@ -2140,7 +2271,7 @@ prepareData2fit <- function(countMatrix, metadata, normalization = TRUE , respon
   }
 
   dtf_countsLong <- countMatrix_2longDtf(countMatrix, response_name)
-  metadata_columnForjoining <- .getColumnWithSampleID(dtf_countsLong, metadata)
+  metadata_columnForjoining <- getColumnWithSampleID(dtf_countsLong, metadata)
   if (is.na(metadata_columnForjoining)) {
     stop("SampleIDs do not seem to correspond between countMatrix and metadata")
   }
@@ -2226,7 +2357,7 @@ test_that("getColumnWithSampleID returns column name with sampleID", {
   expected_output <- "sampleID"
   
   # Get column name with sampleID
-  column_name <- .getColumnWithSampleID(dtf_countLong, mock_data$metadata)
+  column_name <- getColumnWithSampleID(dtf_countLong, mock_data$metadata)
   
   # Check if the output matches the expected output
   expect_identical(column_name, expected_output)
@@ -2280,7 +2411,7 @@ test_that("Throws an error when all-zero genes are encountered", {
 
 ```
 
-```{r functionFitModel, filename = "fitModel"}
+```{r function-fitmodel, filename = "fitmodel"}
 #' Check if Data is Valid for Model Fitting
 #'
 #' This function checks whether the provided data contains all the variables required in the model formula for fitting.
@@ -2294,7 +2425,6 @@ test_that("Throws an error when all-zero genes are encountered", {
 #' data(iris)
 #' formula <- Sepal.Length ~ Sepal.Width + Petal.Length
 #' isValidInput2fit(iris, formula) # Returns TRUE if all required variables are present
-#' @keywords internal
 #' @export
 isValidInput2fit <- function(data2fit, formula){
   vec_bool <- all.vars(formula) %in% colnames(data2fit)
@@ -2324,7 +2454,7 @@ isValidInput2fit <- function(data2fit, formula){
 #' # Drop the random effects related to 'group'
 #' modified_formula <- drop_randfx(formula)
 #'
-#' @importFrom stats terms
+#' @importFrom stats terms drop.terms
 #' @importFrom stats update
 #'
 #' @export
@@ -2394,8 +2524,8 @@ is_fullrank <- function(metadata, formula) {
 #' @return Fitted model object or NULL if there was an error
 #' @export
 #' @examples
-#' .fitModel(formula = mpg ~ cyl + disp, data = mtcars)
-.fitModel <- function(formula, data, ...) {
+#' fitModel(formula = mpg ~ cyl + disp, data = mtcars)
+fitModel <- function(formula, data, ...) {
   # Fit the model using glm.nb from the GLmmTMB package
   model <- glmmTMB::glmmTMB(formula, ..., data = data ) 
   model$frame <- data
@@ -2421,12 +2551,12 @@ is_fullrank <- function(metadata, formula) {
 #' @return Fitted model object or NULL if there was an error
 #' @export
 #' @examples
-#' .subsetData_andfit(group = "setosa", group_by = "Species", 
+#' subsetData_andfit(group = "setosa", group_by = "Species", 
 #'                  formula = Sepal.Length ~ Sepal.Width + Petal.Length, 
 #'                  data = iris )
-.subsetData_andfit <- function(group, group_by, formula, data, ...) {
+subsetData_andfit <- function(group, group_by, formula, data, ...) {
   subset_data <- data[data[[group_by]] == group, ]
-  fit_res <- .fitModel(formula, subset_data, ...)
+  fit_res <- fitModel(formula, subset_data, ...)
   #glance_df <- glance.negbin(group_by ,group , fit_res)
   #tidy_df <- tidy.negbin(group_by ,group,fit_res )
   #list(glance = glance_df, summary = tidy_df)
@@ -2455,7 +2585,7 @@ launchFit <- function(group, group_by, formula, data, ...) {
   tryCatch(
     expr = {
       withCallingHandlers(
-          .subsetData_andfit(group, group_by, formula, data, ...),
+          subsetData_andfit(group, group_by, formula, data, ...),
           warning = function(w) {
             message(paste(Sys.time(), "warning for group", group, ":", conditionMessage(w)))
             invokeRestart("muffleWarning")
@@ -2485,14 +2615,14 @@ launchFit <- function(group, group_by, formula, data, ...) {
 #' @importFrom stats setNames
 #' @export
 #' @examples
-#' .parallel_fit(group_by = "Species", "setosa", 
+#' parallel_fit(group_by = "Species", "setosa", 
 #'                formula = Sepal.Length ~ Sepal.Width + Petal.Length, 
 #'                data = iris, n.cores = 1, log_file = "log.txt" )
-.parallel_fit <- function(groups, group_by, formula, data, n.cores = NULL, log_file,  ...) {
+parallel_fit <- function(groups, group_by, formula, data, n.cores = NULL, log_file,  ...) {
   if (is.null(n.cores)) n.cores <- parallel::detectCores()
   
   clust <- parallel::makeCluster(n.cores, outfile = log_file)
-  parallel::clusterExport(clust, c(".subsetData_andfit", ".fitModel"),  envir=environment())
+  parallel::clusterExport(clust, c("subsetData_andfit", "fitModel"),  envir=environment())
   results_fit <- parallel::parLapply(clust, X = stats::setNames(groups, groups), fun = launchFit, 
                       group_by = group_by, formula = formula, data = data, ...)
                                      
@@ -2525,7 +2655,7 @@ fitModelParallel <- function(formula, data, group_by, n.cores = NULL, log_file =
   
   groups <- unique(data[[group_by]])
   # Fit models in parallel and capture the results
-  results <- .parallel_fit(groups, group_by, formula, data, n.cores, log_file, ...)
+  results <- parallel_fit(groups, group_by, formula, data, n.cores, log_file, ...)
   #results <- mergeListDataframes(results)
   return(results)
 }
@@ -2534,7 +2664,7 @@ fitModelParallel <- function(formula, data, group_by, n.cores = NULL, log_file =
 ```
 
 
-```{r  test-fitData}
+```{r  test-fitmodel}
 
 
 test_that("isValidInput2fit returns TRUE for valid data", {
@@ -2551,28 +2681,28 @@ test_that("isValidInput2fit raises an error for missing variable", {
   expect_error(isValidInput2fit(iris, formula), "Variable NonExistentVariable not found")
 })
 
-test_that(".fitModel returns a fitted model object", {
+test_that("fitModel returns a fitted model object", {
   data(mtcars)
   formula <- mpg ~ cyl + disp
-  fitted_model <- suppressWarnings(.fitModel(formula, mtcars))
-  #expect_warning(.fitModel(formula, mtcars))
+  fitted_model <- suppressWarnings(fitModel(formula, mtcars))
+  #expect_warning(fitModel(formula, mtcars))
   expect_s3_class(fitted_model, "glmmTMB")
   
   # Test with invalid formula
   invalid_formula <- mpg ~ cyl + disp + invalid_var
-  expect_error(.fitModel(invalid_formula, mtcars))
+  expect_error(fitModel(invalid_formula, mtcars))
   
   
    # Additional parameters: 
    #change family + formula
   formula <- Sepal.Length ~ Sepal.Width + Petal.Length + (1 | Species)
-  fitted_models <- suppressWarnings(.fitModel(formula = formula, 
+  fitted_models <- suppressWarnings(fitModel(formula = formula, 
                                                     data = iris, 
                                                     family = glmmTMB::nbinom1(link = "log") ))
   expect_s3_class(fitted_models$call$family, "family")
   expect_equal(fitted_models$call$formula, formula)
   #change control settings
-  fitted_models <- suppressWarnings(.fitModel(formula = formula, 
+  fitted_models <- suppressWarnings(fitModel(formula = formula, 
                                                     data = iris, 
                                                     family = glmmTMB::nbinom1(link = "log"), 
                                                 control = glmmTMB::glmmTMBControl(optCtrl=list(iter.max=1e3,
@@ -2657,23 +2787,23 @@ test_that("Identify full-rank model matrix (with random eff)", {
 #  expect_null(fitted_model_error)
 #})
 
-test_that(".subsetData_andfit returns a glmTMB obj", {
+test_that("subsetData_andfit returns a glmTMB obj", {
   data(iris)
   group <- "setosa"
   group_by <- "Species"
   formula <- Sepal.Length ~ Sepal.Width + Petal.Length
-  fitted_model <- .subsetData_andfit(group, group_by, formula, iris)
+  fitted_model <- subsetData_andfit(group, group_by, formula, iris)
   expect_s3_class(fitted_model, "glmmTMB")
 
   # Test with invalid formula
   invalid_formula <- Sepal.Length ~ Sepal.Width + Petal.Length +  invalid_var
-  expect_error(.subsetData_andfit(group, group_by, invalid_formula, mtcars))
+  expect_error(subsetData_andfit(group, group_by, invalid_formula, mtcars))
   
   
     # Additional parameters: 
    #change family + formula
   formula <- Sepal.Length ~ Sepal.Width + Petal.Length + (1 | Species)
-  fitted_models <- suppressWarnings(.subsetData_andfit(group,
+  fitted_models <- suppressWarnings(subsetData_andfit(group,
                                                        group_by,
                                                        formula = formula, 
                                                         data = iris, 
@@ -2681,7 +2811,7 @@ test_that(".subsetData_andfit returns a glmTMB obj", {
   expect_s3_class(fitted_models$call$family, "family")
   expect_equal(fitted_models$call$formula, formula)
   #change control settings
-  fitted_models <- suppressWarnings(.subsetData_andfit(group,
+  fitted_models <- suppressWarnings(subsetData_andfit(group,
                                                        group_by,
                                                        formula = formula, 
                                                         data = iris, 
@@ -2727,18 +2857,18 @@ test_that("launchFit handles warnings and errors during the fitting process", {
   expect_equal(fitted_models$call$control,  glmmTMB::glmmTMBControl(optCtrl=list(iter.max=1e3,eval.max=1e3)))
 })
 
-test_that(".parallel_fit returns a list of fitted model objects or NULL for any errors", {
+test_that("parallel_fit returns a list of fitted model objects or NULL for any errors", {
   data(iris)
   groups <- unique(iris$Species)
   group_by <- "Species"
   formula <- Sepal.Length ~ Sepal.Width + Petal.Length
-  fitted_models <- .parallel_fit(groups, group_by, formula, iris, log_file = "log.txt", n.cores = 1)
+  fitted_models <- parallel_fit(groups, group_by, formula, iris, log_file = "log.txt", n.cores = 1)
   expect_s3_class(fitted_models$setosa, "glmmTMB")
   expect_length(fitted_models, length(groups))
 
   # Test with invalid formula
   invalid_formula <- blabla ~ cyl + disp 
-  result <- suppressWarnings(.parallel_fit(groups, group_by, invalid_formula,  
+  result <- suppressWarnings(parallel_fit(groups, group_by, invalid_formula,  
                                            iris, log_file = "log.txt",  n.cores = 1))
   expect_equal(result, expected = list(setosa = NULL, versicolor = NULL, virginica = NULL))
   
@@ -2746,7 +2876,7 @@ test_that(".parallel_fit returns a list of fitted model objects or NULL for any
    # Additional parameters: 
    #change family + formula
   formula <- Sepal.Length ~ Sepal.Width + Petal.Length
-  fitted_models <- suppressWarnings(.parallel_fit(formula = formula, 
+  fitted_models <- suppressWarnings(parallel_fit(formula = formula, 
                                                     data = iris, 
                                                     group_by = group_by, 
                                                     groups = "setosa",
@@ -2756,7 +2886,7 @@ test_that(".parallel_fit returns a list of fitted model objects or NULL for any
   expect_s3_class(fitted_models$setosa$call$family, "family")
   expect_equal(fitted_models$setosa$call$formula, formula)
   #change control settings
-  fitted_models <- suppressWarnings(.parallel_fit(formula = formula, 
+  fitted_models <- suppressWarnings(parallel_fit(formula = formula, 
                                                     data = iris, 
                                                     group_by = group_by, 
                                                     groups = "setosa",
@@ -2807,15 +2937,15 @@ test_that("fitModelParallel fits models in parallel for each group and returns a
 ```
 
 
-```{r functionUpdateFitModel, filename = "updateFitModel"}
+```{r function-update_fittedmodel, filename = "update_fittedmodel"}
 
 
-#' Update GLMNB models in parallel.
+#' Update glmmTMB models in parallel.
 #'
-#' This function fits GLMNB models in parallel using multiple cores, allowing for faster computation.
+#' This function fits glmmTMB models in parallel using multiple cores, allowing for faster computation.
 #'
 #' @param formula Formula for the GLMNB model.
-#' @param l_tmb List of GLMNB objects.
+#' @param list_tmb List of glmmTMB objects.
 #' @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.
 #' @param ... Additional arguments to be passed to the glmmTMB::glmmTMB function.
@@ -2830,24 +2960,24 @@ 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, l_tmb, n.cores = NULL, log_file = "log.txt", ...) {
+updateParallel <- function(formula, list_tmb, n.cores = NULL, log_file = "log.txt", ...) {
     
-    isValidInput2fit(l_tmb[[1]]$frame, formula)
+    isValidInput2fit(list_tmb[[1]]$frame, formula)
   
-    is_fullrank(l_tmb[[1]]$frame, formula)
+    is_fullrank(list_tmb[[1]]$frame, formula)
     
     # Fit models update in parallel and capture the results
-    results <- .parallel_update(formula, l_tmb, n.cores, log_file, ...)
+    results <- parallel_update(formula, list_tmb, n.cores, log_file, ...)
     return(results)
 }
 
 
-#' Internal function to fit GLMNB models in parallel.
+#' Internal function to fit glmmTMB models in parallel.
 #'
-#' This function is used internally by \code{\link{updateParallel}} to fit GLMNB models in parallel.
+#' This function is used internally by \code{\link{updateParallel}} to fit glmmTMB models in parallel.
 #'
 #' @param formula Formula for the GLMNB model.
-#' @param l_tmb List of GLMNB objects.
+#' @param list_tmb List of glmmTMB objects.
 #' @param n.cores Number of cores to use for parallel processing.
 #' @param log_file File path for the log output.
 #' @param ... Additional arguments to be passed to the glmmTMB::glmmTMB function.
@@ -2860,13 +2990,13 @@ updateParallel <- function(formula, l_tmb, n.cores = NULL, log_file = "log.txt",
 #' formula <- Sepal.Length ~ Sepal.Width + Petal.Length
 #' fitted_models <- fitModelParallel(formula, iris, group_by, n.cores = 1)
 #' new_formula <- Sepal.Length ~ Sepal.Width 
-#' results <- .parallel_update(new_formula, fitted_models, n.cores = 1)
-.parallel_update <- function(formula, l_tmb, n.cores = NULL, log_file = "log.txt",  ...) {
+#' results <- parallel_update(new_formula, fitted_models, n.cores = 1)
+parallel_update <- function(formula, list_tmb, n.cores = NULL, log_file = "log.txt",  ...) {
   if (is.null(n.cores)) n.cores <- parallel::detectCores()
   clust <- parallel::makeCluster(n.cores, outfile = log_file)
   #l_geneID <- attributes(l_tmb)$names
   parallel::clusterExport(clust, c("launchUpdate", "fitUpdate"),  envir=environment())
-  updated_res <- parallel::parLapply(clust, X = l_tmb, fun = launchUpdate , formula = formula, ...)
+  updated_res <- parallel::parLapply(clust, X = list_tmb, fun = launchUpdate , formula = formula, ...)
   parallel::stopCluster(clust)
   #closeAllConnections()
   return(updated_res)
@@ -2877,7 +3007,7 @@ updateParallel <- function(formula, l_tmb, n.cores = NULL, log_file = "log.txt",
 #'
 #' This function fits and updates a GLMNB model using the provided formula.
 #'
-#' @param glmnb_obj A GLMNB object to be updated.
+#' @param glmm_obj A glmmTMB object to be updated.
 #' @param formula Formula for the updated GLMNB model.
 #' @param ... Additional arguments to be passed to the glmmTMB::glmmTMB function.
 #' @export
@@ -2891,9 +3021,9 @@ updateParallel <- function(formula, l_tmb, n.cores = NULL, log_file = "log.txt",
 #' fitted_models <- fitModelParallel(formula, iris, group_by, n.cores = 1)
 #' new_formula <- Sepal.Length ~ Sepal.Width 
 #' updated_model <- fitUpdate(fitted_models[[1]], new_formula)
-fitUpdate <- function(glmnb_obj, formula , ...){
-  data = glmnb_obj$frame
-  resUpdt <- stats::update(glmnb_obj, formula, ...)
+fitUpdate <- function(glmm_obj, formula , ...){
+  data = glmm_obj$frame
+  resUpdt <- stats::update(glmm_obj, formula, ...)
   resUpdt$frame <- data
   ## family in ... => avoid error in future update
   additional_args <- list(...)
@@ -2910,7 +3040,7 @@ fitUpdate <- function(glmnb_obj, formula , ...){
 #'
 #' This function launches the update process for a GLMNB model, capturing and handling warnings and errors.
 #'
-#' @param glmnb_obj A GLMNB object to be updated.
+#' @param glmm_obj A glmmTMB object to be updated.
 #' @param formula Formula for the updated GLMNB model.
 #' @param ... Additional arguments to be passed to the glmmTMB::glmmTMB function.
 #' @export
@@ -2924,12 +3054,12 @@ fitUpdate <- function(glmnb_obj, formula , ...){
 #' fitted_models <- fitModelParallel(formula, iris, group_by, n.cores = 1)
 #' new_formula <- Sepal.Length ~ Sepal.Width 
 #' updated_model <- launchUpdate(fitted_models[[1]], new_formula)
-launchUpdate <- function(glmnb_obj, formula,  ...) {
-  group = deparse(substitute(glmnb_obj))
+launchUpdate <- function(glmm_obj, formula,  ...) {
+  group = deparse(substitute(glmm_obj))
   tryCatch(
     expr = {
       withCallingHandlers(
-        fitUpdate(glmnb_obj, formula, ...),
+        fitUpdate(glmm_obj, formula, ...),
         warning = function(w) {
           message(paste(Sys.time(), "warning for group", group ,":", conditionMessage(w)))
           invokeRestart("muffleWarning")
@@ -2945,7 +3075,7 @@ launchUpdate <- function(glmnb_obj, formula,  ...) {
 ```
 
 
-```{r  test-updateFit}
+```{r  test-update_fittedmodel}
 # Test updateParallel function
 test_that("updateParallel function returns correct results", {
   # Load the required data
@@ -2967,14 +3097,14 @@ test_that("updateParallel function returns correct results", {
   # Additional parameters: 
    #change family + formula
   new_formula <- Sepal.Length ~ Sepal.Width 
-  updated_model <- suppressWarnings(updateParallel(l_tmb = fitted_models, 
+  updated_model <- suppressWarnings(updateParallel(fitted_models, 
                                                     formula = new_formula,
                                                     n.cores = 1,
                                                     family = glmmTMB::nbinom1(link = "log") ))
   expect_s3_class(updated_model$setosa$call$family, "family")
   expect_equal(updated_model$setosa$call$formula, new_formula)
   #change control settings
-  updated_model <- suppressWarnings(updateParallel(l_tmb = fitted_models, 
+  updated_model <- suppressWarnings(updateParallel(fitted_models, 
                                                  formula = new_formula, 
                                                  family = glmmTMB::nbinom1(link = "log"), 
                                                   n.cores = 1,
@@ -2983,15 +3113,15 @@ test_that("updateParallel function returns correct results", {
   expect_equal(updated_model$setosa$call$control,  glmmTMB::glmmTMBControl(optCtrl=list(iter.max=1e3,eval.max=1e3)))
   
   # Update an updated model
-  updated_updated_model <- suppressWarnings(updateParallel(l_tmb = updated_model, 
+  updated_updated_model <- suppressWarnings(updateParallel(updated_model, 
                                                  formula = new_formula, 
                                                   n.cores = 1,
                                                  family = glmmTMB::ziGamma(link = "inverse")))
   expect_s3_class(updated_updated_model$setosa$call$family,  "family")
 })
 
-# Test .parallel_update function
-test_that(".parallel_update function returns correct results", {
+# Test parallel_update function
+test_that("parallel_update function returns correct results", {
 # Load the required data
   data(iris)
   groups <- unique(iris$Species)
@@ -2999,14 +3129,14 @@ test_that(".parallel_update function returns correct results", {
   formula <- Sepal.Length ~ Sepal.Width + Petal.Length
   fitted_models <- fitModelParallel(formula, iris, group_by, n.cores = 1)
   new_formula <- Sepal.Length ~ Sepal.Width 
-  results <- .parallel_update(new_formula, fitted_models, n.cores = 1)
+  results <- parallel_update(new_formula, fitted_models, n.cores = 1)
   expect_is(results, "list")
   expect_equal(length(results), length(fitted_models))
   expect_is(results$setosa, "glmmTMB")
 
   #uncorrect formula 
   new_formula <- Sepal.Length ~ blabla
-  results <- .parallel_update(new_formula, fitted_models, n.cores = 1)
+  results <- parallel_update(new_formula, fitted_models, n.cores = 1)
   expect_is(results, "list")
   expect_equal(length(results), length(fitted_models))
   expect_equal(results$setosa, NULL)
@@ -3014,14 +3144,14 @@ test_that(".parallel_update function returns correct results", {
   # Additional parameters: 
    #change family + formula
   new_formula <- Sepal.Length ~ Sepal.Width 
-  updated_model <- suppressWarnings(.parallel_update(l_tmb = fitted_models, 
+  updated_model <- suppressWarnings(parallel_update(fitted_models, 
                                                      formula = new_formula,
                                                       n.cores = 1,
                                                       family = glmmTMB::nbinom1(link = "log") ))
   expect_s3_class(updated_model$setosa$call$family, "family")
   expect_equal(updated_model$setosa$call$formula, new_formula)
   #change control
-  updated_model <- suppressWarnings(.parallel_update(l_tmb = fitted_models, 
+  updated_model <- suppressWarnings(parallel_update(fitted_models, 
                                                  formula = new_formula, 
                                                   n.cores = 1,
                                                  family = glmmTMB::nbinom1(link = "log"), 
@@ -3084,7 +3214,7 @@ test_that("launchUpdate function returns correct results", {
 
 ```
 
-```{r functionTidyGLM, filename = "tidy_glmmTMB"}
+```{r function-tidy_glmmTMB, filename = "tidy_glmmTMB"}
 
 
 #' Extract Fixed Effects from a GLMMTMB Model Summary
@@ -3155,7 +3285,7 @@ getTidyGlmmTMB <- function(glm_TMB, ID){
 #'
 #' This function takes a list of glmmTMB models and extracts a tidy summary of the fixed and random effects from each model. It then combines the results into a single data frame.
 #'
-#' @param l_tmb A list of glmmTMB model objects.
+#' @param list_tmb A list of glmmTMB model objects.
 #' @return A data frame containing a tidy summary of the fixed and random effects from all glmmTMB models in the list.
 #' @export
 #' @examples
@@ -3163,13 +3293,13 @@ getTidyGlmmTMB <- function(glm_TMB, ID){
 #' model2 <- glmmTMB::glmmTMB(Petal.Length ~ Sepal.Length + Sepal.Width + (1 | Species), data = iris)
 #' model_list <- list(Model1 = model1, Model2 = model2)
 #' tidy_summary <- tidy_tmb(model_list)
-tidy_tmb <- function(l_tmb){
-    if (identical(class(l_tmb), "glmmTMB")) return(getTidyGlmmTMB(l_tmb, "glmmTMB"))
-    attributes(l_tmb)$names
-    l_tidyRes <- lapply(attributes(l_tmb)$names,
+tidy_tmb <- function(list_tmb){
+    if (identical(class(list_tmb), "glmmTMB")) return(getTidyGlmmTMB(list_tmb, "glmmTMB"))
+    attributes(list_tmb)$names
+    l_tidyRes <- lapply(attributes(list_tmb)$names,
                  function(ID)
                    {
-                      glm_TMB <- l_tmb[[ID]]
+                      glm_TMB <- list_tmb[[ID]]
                       getTidyGlmmTMB(glm_TMB, ID)
                   }
                 )
@@ -3202,22 +3332,6 @@ build_missingColumn_with_na <- function(df, l_columns = c("effect", "component",
   return(df)
 }
 
-#' Remove Duplicated Words from Strings
-#'
-#' This function takes a vector of strings and removes duplicated words within each string.
-#'
-#' @param strings A character vector containing strings with potential duplicated words.
-#' @return A character vector with duplicated words removed from each string.
-#' @export
-#' @examples
-#'
-#' words <- c("hellohello", "worldworld", "programmingprogramming", "R isis great")
-#' cleaned_words <- removeDuplicatedWord(words)
-removeDuplicatedWord <- function(strings){
-  gsub("(.*)\\1+", "\\1", strings, perl = TRUE)
-}
-
-
 
 
 #' Convert Correlation Matrix to Data Frame
@@ -3358,35 +3472,10 @@ renameColumns <- function(df, old_names  = c("Estimate","Std..Error", "z.value",
   return(df)
 }
 
-
-
-#' Reorder the columns of a dataframe
-#'
-#' This function reorders the columns of a dataframe according to the specified column order.
-#'
-#' @param df The input dataframe.
-#' @param columnOrder A vector specifying the desired order of columns.
-#'
-#' @return A dataframe with columns reordered according to the specified column order.
-#' @export
-#' @examples
-#' # Example dataframe
-#' df <- data.frame(A = 1:3, B = 4:6, C = 7:9)
-#'
-#' # Define the desired column order
-#' columnOrder <- c("B", "C", "A")
-#'
-#' # Reorder the columns of the dataframe
-#' df <- reorderColumns(df, columnOrder)
-reorderColumns <- function(df, columnOrder) {
-  df <- df[, columnOrder, drop = FALSE]
-  return(df)
-}
-
 ```
 
 
-```{r  test-tidyGLM}
+```{r  test-tidy_glmmTMB}
 
 test_that("extract_fixed_effect returns the correct results for glmmTMB models", {
   data(iris)
@@ -3440,13 +3529,6 @@ test_that("build_missingColumn_with_na returns the correct results", {
 })
 
 
-test_that("removeDuplicatedWord returns expected output", {
-  words <- c("hellohello", "worldworld", "programmingprogramming", "R isis great")
-  cleaned_words <- removeDuplicatedWord(words)
-  expect_equal(cleaned_words, c("hello", "world", "programming", "R is great"))
-})
-
-
 
 test_that("correlation_matrix_2df returns expected output",{
 
@@ -3494,19 +3576,6 @@ test_that("renameColumns returns expected output",{
   expect_equal(dim(renamed_df), dim(df))
 })
     
-
-test_that("reorderColumns returns expected output",{
-    df <- data.frame(A = 1:3, B = 4:6, C = 7:9)
-    # Define the desired column order
-    columnOrder <- c("B", "C", "A")
-    # Reorder the columns of the dataframe
-    df_reorder <- reorderColumns(df, columnOrder)
-    expect_equal(colnames(df_reorder), columnOrder)
-    expect_equal(dim(df_reorder), dim(df))
-
-})
-
-
 test_that("tidy_tmb returns expected output",{
   model1 <- glmmTMB::glmmTMB(Sepal.Length ~ Sepal.Width + Petal.Length + (1 | Species), data = iris)
   model2 <- glmmTMB::glmmTMB(Petal.Length ~ Sepal.Length + Sepal.Width + (1 | Species), data = iris)
@@ -3543,15 +3612,14 @@ test_that("tidy_tmb returns expected output",{
 })
 ```
 
-
-```{r functionGlanceGLM, filename = "glance_tmb"}
+```{r function-glance_glmmTMB, filename = "glance_glmmTMB"}
 
 #' Extracts the summary statistics from a list of glmmTMB models.
 #'
 #' This function takes a list of glmmTMB models and extracts the summary statistics (AIC, BIC, logLik, deviance,
 #' df.resid, and dispersion) for each model and returns them as a single DataFrame.
 #'
-#' @param l_tmb A list of glmmTMB models or a unique glmmTMB obj model
+#' @param list_tmb A list of glmmTMB models or a unique glmmTMB obj model
 #' @return A DataFrame with the summary statistics for all the glmmTMB models in the list.
 #' @export
 #' @importFrom stats setNames
@@ -3560,10 +3628,10 @@ test_that("tidy_tmb returns expected output",{
 #' models <-  fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length, 
 #'                            group_by = "Species",n.cores = 1, data = iris)
 #' result <- glance_tmb(models)
-glance_tmb <- function(l_tmb){
-  if (identical(class(l_tmb), "glmmTMB")) return(getGlance(l_tmb))
-  l_group <- attributes(l_tmb)$names
-  l_glance <- lapply(stats::setNames(l_group, l_group), function(group) getGlance(l_tmb[[group]]))
+glance_tmb <- function(list_tmb){
+  if (identical(class(list_tmb), "glmmTMB")) return(getGlance(list_tmb))
+  l_group <- attributes(list_tmb)$names
+  l_glance <- lapply(stats::setNames(l_group, l_group), function(group) getGlance(list_tmb[[group]]))
   return(do.call("rbind", l_glance))
 }
 
@@ -3592,7 +3660,7 @@ getGlance <- function(x){
 ```
 
 
-```{r testGlanceGLM }
+```{r test-glance_glmmTMB }
 
 test_that("glance_tmb returns the summary statistics for multiple models", {
   data(iris)
@@ -3631,7 +3699,7 @@ test_that("getGlance returns the summary statistics for a single model", {
 ```
 
 
-```{r functionPlotMetrics, filename = "plot_metrics"}
+```{r function-plot_metrics, filename = "plot_metrics"}
 
 #' Subset the glance DataFrame based on selected variables.
 #'
@@ -3667,7 +3735,7 @@ subset_glance <- function(glance_df, focus){
 #' This function generates a density plot of the specified metrics for the given
 #' list of generalized linear mixed models (GLMMs).
 #'
-#' @param l_tmb A list of GLMM objects to extract metrics from.
+#' @param list_tmb A list of GLMM objects to extract metrics from.
 #' @param focus A character vector specifying the metrics to focus on. Possible
 #'   values include "AIC", "BIC", "logLik", "deviance", "df.resid", and
 #'   "dispersion". If \code{NULL}, all available metrics will be plotted.
@@ -3683,8 +3751,8 @@ subset_glance <- function(glance_df, focus){
 #' models_list <-  fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length, 
 #'                      group_by = "Species",n.cores = 1, data = iris)
 #' metrics_plot(models_list, focus = c("AIC", "BIC", "deviance"))
-metrics_plot <- function(l_tmb, focus = NULL) {
-  glance_df <- glance_tmb(l_tmb)
+metrics_plot <- function(list_tmb, focus = NULL) {
+  glance_df <- glance_tmb(list_tmb)
   glance_df$group_id <- rownames(glance_df)
   if (!is.null(focus)) {
     glance_df <- subset_glance(glance_df, focus)
@@ -3702,7 +3770,7 @@ metrics_plot <- function(l_tmb, focus = NULL) {
 
 ```
 
-```{r testPlotMetrics }
+```{r test-plot_metrics }
 
 
 test_that("subset_glance subsets the glance DataFrame correctly", {
@@ -3745,7 +3813,7 @@ test_that("metrics_plot returns a ggplot object", {
 
 
 
-```{r functionEvalDispersion, filename = "evaluateDispersion"}
+```{r function-evaluate_dispersion, filename = "evaluate_dispersion"}
 
 #' Evaluate Dispersion Comparison
 #'
@@ -3786,7 +3854,6 @@ evaluateDispersion <- function(TMB_dispersion_df, DESEQ_dispersion_df, color2use
 #' @examples
 #' \dontrun{
 #' dispersion_comparison <- getDispersionComparison(inferred_disp, actual_disp)
-#' print(dispersion_comparison)
 #' }
 getDispersionComparison <- function(inferred_dispersion, actual_dispersion) {
   actual_disp <- data.frame(actual_dispersion = actual_dispersion)
@@ -3801,7 +3868,7 @@ getDispersionComparison <- function(inferred_dispersion, actual_dispersion) {
 #'
 #' Extracts inferred dispersion values from a DESeq2 wrapped object.
 #'
-#' @param deseq_wrapped A DESeq2 wrapped object containing dispersion values.
+#' @param dds_wrapped A DESeq2 wrapped object containing dispersion values.
 #'
 #' @return A data frame containing inferred dispersion values.
 #' 
@@ -3809,11 +3876,10 @@ getDispersionComparison <- function(inferred_dispersion, actual_dispersion) {
 #'
 #' @examples
 #' \dontrun{
-#' dispersion_df <- extractDESeqDispersion(deseq2_object)
-#' print(dispersion_df)
+#' dispersion_df <- extract_ddsDispersion(deseq2_object)
 #' }
-extractDESeqDispersion <- function(deseq_wrapped) {
-  inferred_dispersion <- data.frame(inferred_dispersion = deseq_wrapped$dispersion)
+extract_ddsDispersion <- function(dds_wrapped) {
+  inferred_dispersion <- data.frame(inferred_dispersion = dds_wrapped$dispersion)
   inferred_dispersion$geneID <- rownames(inferred_dispersion)
   rownames(inferred_dispersion) <- NULL
   return(inferred_dispersion)
@@ -3824,7 +3890,7 @@ extractDESeqDispersion <- function(deseq_wrapped) {
 #'
 #' Extracts inferred dispersion values from a TMB result object.
 #'
-#' @param l_tmb A TMB result object containing dispersion values.
+#' @param list_tmb A TMB result object containing dispersion values.
 #'
 #' @return A data frame containing inferred dispersion values.
 #' 
@@ -3832,11 +3898,10 @@ extractDESeqDispersion <- function(deseq_wrapped) {
 #'
 #' @examples
 #' \dontrun{
-#' dispersion_df <- extractTMBDispersion(tmb_result)
-#' print(dispersion_df)
+#' dispersion_df <- extract_tmbDispersion(tmb_result)
 #' }
-extractTMBDispersion <- function(l_tmb) {
-  glanceRes <- glance_tmb(l_tmb)
+extract_tmbDispersion <- function(list_tmb) {
+  glanceRes <- glance_tmb(list_tmb)
   inferred_dispersion <- data.frame(inferred_dispersion = glanceRes$dispersion)
   inferred_dispersion$geneID <- rownames(glanceRes)
   rownames(inferred_dispersion) <- NULL
@@ -3880,7 +3945,7 @@ dispersion_plot <- function(eval_dispersion, ...) {
 
 ```
 
-```{r testPlotMetrics }
+```{r test-evaluate_dispersion }
 
 
 # Example data
@@ -3897,7 +3962,7 @@ test_that("dispersion_plot function works correctly", {
   expect_s3_class(disp_plot, "gg")
 })
 
-test_that("extractTMBDispersion function extracts dispersion correctly", {
+test_that("extract_tmbDispersion function extracts dispersion correctly", {
    N_GENES = 100
   MAX_REPLICATES = 5
   MIN_REPLICATES = 5
@@ -3908,11 +3973,11 @@ test_that("extractTMBDispersion function extracts dispersion correctly", {
   l_res <- fitModelParallel(formula = kij ~ varA,
                           data = data2fit, group_by = "geneID",
                           family = glmmTMB::nbinom2(link = "log"), n.cores = 1)
-  extracted_disp <- extractTMBDispersion(l_res)
+  extracted_disp <- extract_tmbDispersion(l_res)
   expect_identical(colnames(extracted_disp), c("inferred_dispersion", "geneID"))
 })
 
-test_that("extractDESeqDispersion function extracts dispersion correctly", {
+test_that("extract_ddsDispersion function extracts dispersion correctly", {
    N_GENES = 100
   MAX_REPLICATES = 5
   MIN_REPLICATES = 5
@@ -3924,9 +3989,9 @@ test_that("extractDESeqDispersion function extracts dispersion correctly", {
       colData = mock_data$metadata,
       design = ~ varA)
   dds <- DESeq2::DESeq(dds, quiet = TRUE)
-  deseq_wrapped = wrapper_DESeq2(dds, 2, "greaterAbs")
+  deseq_wrapped = wrap_dds(dds, 2, "greaterAbs")
   
-  extracted_disp <- extractDESeqDispersion(deseq_wrapped)
+  extracted_disp <- extract_ddsDispersion(deseq_wrapped)
   expect_identical(colnames(extracted_disp), c("inferred_dispersion", "geneID"))
 })
 
@@ -3942,7 +4007,7 @@ test_that("getDispersionComparison function works correctly", {
                           data = data2fit, group_by = "geneID",
                           family = glmmTMB::nbinom2(link = "log"), n.cores = 1)
   
-  tmb_disp_inferred <- extractTMBDispersion(l_res)
+  tmb_disp_inferred <- extract_tmbDispersion(l_res)
     
   comparison <- getDispersionComparison(tmb_disp_inferred, mock_data$groundTruth$gene_dispersion)
   expect_identical(colnames(comparison), c("actual_dispersion",  "geneID", "inferred_dispersion"))
@@ -3964,12 +4029,12 @@ test_that("evaluateDispersion function works correctly", {
       colData = mock_data$metadata,
       design = ~ varA)
   dds <- DESeq2::DESeq(dds, quiet = TRUE)
-  deseq_wrapped = wrapper_DESeq2(dds, 2, "greaterAbs")
+  deseq_wrapped = wrap_dds(dds, 2, "greaterAbs")
   
-  tmb_disp_inferred <- extractTMBDispersion(l_res)
+  tmb_disp_inferred <- extract_tmbDispersion(l_res)
   TMB_dispersion_df <- getDispersionComparison(tmb_disp_inferred, mock_data$groundTruth$gene_dispersion)
   TMB_dispersion_df$from <- 'HTRfit'
-  DESEQ_disp_inferred <- extractDESeqDispersion(deseq_wrapped)
+  DESEQ_disp_inferred <- extract_ddsDispersion(deseq_wrapped)
   DESEQ_dispersion_df <- getDispersionComparison(DESEQ_disp_inferred , mock_data$groundTruth$gene_dispersion)
   DESEQ_dispersion_df$from <- 'DESeq2'
     
@@ -3985,7 +4050,7 @@ test_that("evaluateDispersion function works correctly", {
 
 
 
-```{r function-seqDepth, filename =  "scalingSequencingDepth"}
+```{r function-sequencing_depth_scaling, filename =  "sequencing_depth_scaling"}
 
 #' Scale Counts Table
 #'
@@ -4018,7 +4083,7 @@ scaleCountsTable <- function(countsTable, seq_depth){
 
 ```
 
-```{r  test-scalingSequencingDepth}
+```{r  test-sequencing_depth_scaling}
 
 # Test case 1: Scaling with valid min_seq_depth and max_seq_depth
 test_that("Valid scaling of counts table", {
@@ -4040,7 +4105,7 @@ test_that("Valid scaling of counts table", {
 
 
 
-```{r function-geneExpressionScaling, filename =  "scalingGeneExpression"}
+```{r function-basal_expression_scaling, filename =  "basal_expression_scaling"}
 
 
 
@@ -4071,7 +4136,7 @@ getBinExpression <- function(dtf_coef, n_bins){
 
 #' Generate BE data.
 #' 
-#' This function generates BE data for a given number of genes, in a vector of BE values.
+#' This function generates basal expression data for a given number of genes, in a vector of basal expression values.
 #' 
 #' @param n_genes The number of genes to generate BE data for.
 #' @param basal_expression a numeric vector from which sample BE for eacg genes
@@ -4079,10 +4144,10 @@ getBinExpression <- function(dtf_coef, n_bins){
 #' @return A data frame containing gene IDs, BE values
 #' 
 #' @examples
-#' generate_BE(n_genes = 100, 10)
+#' generate_basal_expression(n_genes = 100, 10)
 #' 
 #' @export
-generate_BE <- function(n_genes, basal_expression) {
+generate_basal_expression <- function(n_genes, basal_expression) {
   ## --avoid bug if one value in basal_expr
   pool2sample <- c(basal_expression, basal_expression)
   BE <- sample(x = pool2sample, size = n_genes, replace = T)
@@ -4097,7 +4162,7 @@ generate_BE <- function(n_genes, basal_expression) {
 #'
 #' This function takes the coefficients data frame \code{dtf_coef} and computes
 #' basal expression for gene expression. The scaling factors are generated 
-#' using the function \code{generate_BE}.
+#' using the function \code{generate_basal_expression}.
 #'
 #' @param dtf_coef A data frame containing the coefficients for gene expression.
 #' @param n_genes number of genes in simulation
@@ -4113,7 +4178,7 @@ generate_BE <- function(n_genes, basal_expression) {
 #' dtf_coef <- getLog_qij(dtf_coef)
 #' addBasalExpression(dtf_coef, N_GENES, 1)
 addBasalExpression <- function(dtf_coef, n_genes, basal_expression){
-    BE_df  <-  generate_BE(n_genes, basal_expression )
+    BE_df  <-  generate_basal_expression(n_genes, basal_expression )
     dtf_coef <- join_dtf(dtf_coef, BE_df, "geneID", "geneID")
     return(dtf_coef) 
 }
@@ -4123,17 +4188,17 @@ addBasalExpression <- function(dtf_coef, n_genes, basal_expression){
 
 ```
 
-```{r  test-geneExpressionScaling}
+```{r  test-basal_expression_scaling}
 
-test_that("generate_BE returns correct number of genes", {
-  be_data <- generate_BE(n_genes = 100, 1)
+test_that("generate_basal_expression returns correct number of genes", {
+  be_data <- generate_basal_expression(n_genes = 100, 1)
   expect_equal(nrow(be_data), 100)
 })
 
 
-test_that("generate_BE returns BE values within specified vector", {
+test_that("generate_basal_expression returns BE values within specified vector", {
   BE_vec <- c(1, 2, 33, 0.4)
-  be_data <- generate_BE(n_genes = 100, BE_vec)
+  be_data <- generate_basal_expression(n_genes = 100, BE_vec)
   expect_true(all(be_data$basalExpr %in% BE_vec))
 })
 
@@ -4201,7 +4266,7 @@ test_that("getBinExpression handles negative values correctly", {
 
 
 
-```{r functionActualMainFixEff, filename =  "actualMainFixEffects" }
+```{r function-actual_mainfixeffects, filename =  "actual_mainfixeffects" }
 
 #' Calculate average values by group
 #'
@@ -4221,22 +4286,6 @@ averageByGroup <- function(data, column, group_by) {
   return(result)
 }
 
-#' Convert specified columns to factor
-#'
-#' @param data The input data frame
-#' @param columns The column names to be converted to factors
-#' @return The modified data frame with specified columns converted to factors
-#' @export
-convert2Factor <- function(data, columns) {
-  if (is.character(columns)) {
-    columns <- match(columns, colnames(data))
-  }
-
-  if (length(columns) > 1) data[, columns] <- lapply(data[, columns], as.factor )
-  else data[, columns] <- as.factor(data[, columns])
-  return(data)
-}
-
 #' Subset Fixed Effect Inferred Terms
 #'
 #' This function subsets the tidy TMB object to extract the fixed effect inferred terms
@@ -4442,13 +4491,9 @@ getActualMainFixEff <- function( l_term , fixeEff_dataActual , df_actualIntercep
   return(df_actual)
 }
 
-
-
-
-
 ```
 
-```{r test-actualMainFixEff}
+```{r test-actual_mainfixeffects}
 
 test_that("Test for subsetFixEffectInferred function", {
   # Prepare the test data
@@ -4496,23 +4541,7 @@ test_that("averageByGroup returns correct average values", {
   expect_equal(result$logQij_mean, c(1, 3,5, 7, 2, 4, 6, 8))  # Average values
 })
 
-# Tests for convert2Factor
-test_that("convert2Factor converts specified columns to factors", {
-  # Create a sample data frame
-  data <- data.frame(
-    Category1 = c("A", "B", "A", "B"),
-    Category2 = c("X", "Y", "X", "Z"),
-    Value = 1:4,
-    stringsAsFactors = FALSE
-  )
-  
-  # Convert columns to factors
-  result <- convert2Factor(data, columns = c("Category1", "Category2"))
-  
-  # Check the output
-  expect_is(result$Category1, "factor")  # Category1 column converted to factor
-  expect_is(result$Category2, "factor")  # Category2 column converted to factor
-})
+
 
 # Tests for findAttribute
 test_that("findAttribute returns the correct attribute", {
@@ -4669,10 +4698,9 @@ test_that("generateActualForMainFixEff returns correct values for main fixed eff
   expect_equal(df_term, expected)
 })
 
-
 ```
 
-```{r functionActualInteractionFixEff, filename =  "actualInteractionFixEffects" }
+```{r function-actual_interactionfixeffects, filename =  "actual_interactionfixeffects" }
 #' Filter DataFrame
 #'
 #' Filter a DataFrame based on the specified filter list.
@@ -5003,7 +5031,7 @@ computeActualInteractionFixEff <- function(l_interactionTerm, categorical_vars,
 }
 ```
 
-```{r test-actualInteractionFixEff }
+```{r test-actual_interactionfixeffects }
 
 test_that("filter_dataframe retourne le dataframe filtré correctement", {
   # Créer un exemple de dataframe
@@ -5417,7 +5445,7 @@ inferenceToExpected_withFixedEff <- function(tidy_tmb , df_ground_truth) {
 ```
 
 
-```{r function-waldTest, filename =  "waldTest" }
+```{r function-waldtest, filename =  "waldtest" }
 
 #' Wald test for hypothesis testing
 #'
@@ -5497,7 +5525,7 @@ tidy_results <- function(list_tmb, coeff_threshold = 0, alternative_hypothesis =
 
 
 
-```{r test-waldTest}
+```{r test-waldtest}
 
 # Test unitaires
 test_that("wald_test performs correct tests", {
@@ -5562,7 +5590,7 @@ test_that("results function performs statistical tests correctly", {
 
 
 
-```{r function-rocPlot, filename = "ROCplot"}
+```{r function-roc_plot, filename = "roc_plot"}
 
 
 #' Get Labels for Expected Differential Expression
@@ -5656,7 +5684,7 @@ roc_plot <- function(comparison_df, ...) {
 
 ```
 
-```{r test-rocPlot}
+```{r test-roc_plot}
 
 
 # Test cases for getLabelExpected function
@@ -5710,7 +5738,7 @@ test_that("ROC plot is generated correctly", {
 ```
 
 
-```{r function-countsPlot, filename = "countsPlot"}
+```{r function-counts_plot, filename = "counts_plot"}
 
 #' Generate a density plot of gene counts
 #'
@@ -5743,7 +5771,7 @@ counts_plot <- function(mock_obj){
 
 ```
 
-```{r test-countsPlot}
+```{r test-counts_plot}
 
 
 
@@ -5764,7 +5792,7 @@ test_that("Counts plot is generated correctly", {
 
 
 
-```{r function-identityPlot, filename = "identityPlot"}
+```{r function-identity_plot, filename = "identity_plot"}
 
 #' Generate an identity plot
 #'
@@ -5803,7 +5831,7 @@ identity_plot <- function(comparison_df, ...){
 
 ```
 
-```{r test-identityPlot}
+```{r test-identity_plot}
 
 
 # Test cases
@@ -5824,7 +5852,7 @@ test_that("Identity plot is generated correctly", {
 ```
 
 
-```{r function-simulationReport, filename =  "simulationReport" }
+```{r function-simulation_report, filename =  "simulation_report" }
 
 #' Export the Analysis Report to a File
 #'
@@ -5916,19 +5944,19 @@ simulationReport <- function(mock_obj, list_tmb = NULL,
       TMB_comparison_df <- compareInferenceToExpected(tidyRes, mock_obj$groundTruth$effects, formula_used)
       TMB_comparison_df <- getLabelExpected(TMB_comparison_df, coeff_threshold, alt_hypothesis)
       TMB_comparison_df$from <- "HTRfit"
-      tmb_disp_inferred <- extractTMBDispersion(list_tmb)
+      tmb_disp_inferred <- extract_tmbDispersion(list_tmb)
       TMB_dispersion_df <- getDispersionComparison(tmb_disp_inferred, mock_data$groundTruth$gene_dispersion)
       TMB_dispersion_df$from <- 'HTRfit'
   }
   
   if (!is.null(dds_obj)){
-      deseq2_wrapped <- wrapper_DESeq2(dds, coeff_threshold, alt_hypothesis)
+      deseq2_wrapped <- wrap_dds(dds, coeff_threshold, alt_hypothesis)
       DESEQ_comparison_df <- inferenceToExpected_withFixedEff(deseq2_wrapped$fixEff, mock_obj$groundTruth$effects)
       DESEQ_comparison_df <- getLabelExpected(DESEQ_comparison_df, coeff_threshold, alt_hypothesis)
       DESEQ_comparison_df$from <- "DESeq2"
       DESEQ_comparison_df$component <- NA
       DESEQ_comparison_df$group <- NA
-      DESEQ_disp_inferred <- extractDESeqDispersion(deseq2_wrapped)
+      DESEQ_disp_inferred <- extract_ddsDispersion(deseq2_wrapped)
       DESEQ_dispersion_df <- getDispersionComparison(DESEQ_disp_inferred , mock_data$groundTruth$gene_dispersion)
       DESEQ_dispersion_df$from <- 'DESeq2'
   }
@@ -5963,7 +5991,7 @@ simulationReport <- function(mock_obj, list_tmb = NULL,
 
 
 
-```{r test-simulationReport}
+```{r test-simulation_report}
 
 
 # Test case 1: Testing with a sample data frame
@@ -5993,7 +6021,7 @@ test_that("Handling non-numeric values in the data frame", {
 ```
 
 
-```{r function-deseq2, filename =  "wrapperDESeq2" }
+```{r function-wrapper_dds, filename =  "wrapper_dds" }
 
 #' Wrapper Function for DESeq2 Analysis
 #'
@@ -6023,9 +6051,9 @@ test_that("Handling non-numeric values in the data frame", {
 #' dds <- DESeq2::DESeqDataSetFromMatrix(mock_data$counts , 
 #'                    mock_data$metadata, ~ genotype + environment)
 #' dds <- DESeq2::DESeq(dds, quiet = TRUE)
-#' result <- wrapper_DESeq2(dds, lfcThreshold = 1, altHypothesis = "greater")
+#' result <- wrap_dds(dds, lfcThreshold = 1, altHypothesis = "greater")
 #' @export
-wrapper_DESeq2 <- function(dds, lfcThreshold , altHypothesis, correction_method = "BH") {
+wrap_dds <- function(dds, lfcThreshold , altHypothesis, correction_method = "BH") {
   dds_full <- S4Vectors::mcols(dds) %>% as.data.frame()
   
   ## -- dispersion
@@ -6034,7 +6062,7 @@ wrapper_DESeq2 <- function(dds, lfcThreshold , altHypothesis, correction_method
   names(dispersion) <- rownames(dds_full)
 
   ## -- coeff
-  inference_df <- get_inference(dds_full, lfcThreshold, altHypothesis, correction_method)
+  inference_df <- get_inference_dds(dds_full, lfcThreshold, altHypothesis, correction_method)
   res <- list(dispersion = dispersion, fixEff = inference_df)
   return(res)
 }
@@ -6055,13 +6083,13 @@ wrapper_DESeq2 <- function(dds, lfcThreshold , altHypothesis, correction_method
 #' @examples
 #' \dontrun{
 #' # Example usage of the function
-#' inference_result <- get_inference(dds_full, lfcThreshold = 0.5, 
+#' inference_result <- get_inference_dds(dds_full, lfcThreshold = 0.5, 
 #'                                    altHypothesis = "greater", 
 #'                                    correction_method = "BH")
 #' }
 #' @importFrom stats p.adjust
 #' @export
-get_inference <- function(dds_full, lfcThreshold, altHypothesis, correction_method){
+get_inference_dds <- function(dds_full, lfcThreshold, altHypothesis, correction_method){
 
   ## -- build subdtf
   stdErr_df <- getSE_df(dds_full)
@@ -6154,10 +6182,10 @@ getEstimate_df <- function(dds_full){
 ```
 
 
-```{r test-wrapperDESeq2}
+```{r test-wrapper_dds}
 
 
-test_that("get_inference returns a data frame with correct columns", {
+test_that("get_inference_dds returns a data frame with correct columns", {
   # Create a sample dds_full data frame
   N_GENES = 100
   MAX_REPLICATES = 5
@@ -6172,7 +6200,7 @@ test_that("get_inference returns a data frame with correct columns", {
   dds_full <- S4Vectors::mcols(dds) %>% as.data.frame()
   
   # Call the function
-  inference_results <- get_inference(dds_full, lfcThreshold = 0.5, altHypothesis = "greater", correction_method = "BH")
+  inference_results <- get_inference_dds(dds_full, lfcThreshold = 0.5, altHypothesis = "greater", correction_method = "BH")
   
   # Check if the returned object is a data frame
   expect_true(is.data.frame(inference_results))
@@ -6259,7 +6287,7 @@ test_that("wrapperDESeq2 function works correctly", {
   mock_data <- mock_rnaseq(input_var_list, N_GENES, MIN_REPLICATES, max_replicates = MAX_REPLICATES)
   dds <- DESeq2::DESeqDataSetFromMatrix(mock_data$counts , mock_data$metadata, ~ genotype + environment)
   dds <- DESeq2::DESeq(dds, quiet = TRUE)
-  deseq2_wrapped <- wrapper_DESeq2(dds, 0.2, "greaterAbs")
+  deseq2_wrapped <- wrap_dds(dds, 0.2, "greaterAbs")
   
   expect_true(is.list(deseq2_wrapped))
   
@@ -6283,7 +6311,7 @@ test_that("wrapperDESeq2 function works correctly", {
 #'
 #' This function handles ANOVA errors and warnings during the ANOVA calculation process.
 #'
-#' @param l_TMB A list of fitted glmmTMB models.
+#' @param list_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.
 #' 
@@ -6291,17 +6319,17 @@ test_that("wrapperDESeq2 function works correctly", {
 #' @export
 #' 
 #' @examples
-#' l_tmb <- fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length,
+#' list_tmb <- fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length,
 #'                           data = iris, group_by = "Species", n.cores = 1)
-#' anova_res <- handleAnovaError(l_tmb, "setosa", type = "III")
+#' anova_res <- handleAnovaError(list_tmb, "setosa", type = "III")
 #'
 #' @importFrom car Anova
 #' @export
-handleAnovaError <- function(l_TMB, group, ...) {
+handleAnovaError <- function(list_tmb, group, ...) {
   tryCatch(
     expr = {
       withCallingHandlers(
-        car::Anova(l_TMB[[group]], ...),
+        car::Anova(list_tmb[[group]], ...),
         warning = function(w) {
           message(paste(Sys.time(), "warning for group", group, ":", conditionMessage(w)))
           invokeRestart("muffleWarning")
@@ -6321,7 +6349,7 @@ handleAnovaError <- function(l_TMB, group, ...) {
 #' models in parallel for different groups specified in the list. It returns a list
 #' of ANOVA results for each group.
 #'
-#' @param l_tmb A list of \code{glmmTMB} models, with model names corresponding to the groups.
+#' @param list_tmb 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.
@@ -6329,14 +6357,14 @@ handleAnovaError <- function(l_TMB, group, ...) {
 #' @examples
 #' # Perform ANOVA
 #' data(iris)
-#' l_tmb<- fitModelParallel( Sepal.Length ~ Sepal.Width  + Petal.Length, 
+#' list_tmb<- fitModelParallel( Sepal.Length ~ Sepal.Width  + Petal.Length, 
 #'                          data = iris, group_by = "Species", n.cores = 1 )
-#' anov_res <- anovaParallel(l_tmb , type = "III")
+#' anov_res <- anovaParallel(list_tmb , type = "III")
 #' @importFrom stats anova
 #' @export
-anovaParallel <- function(l_tmb, ...) {
-  l_group <- attributes(l_tmb)$names
-  lapply(stats::setNames(l_group, l_group), function(group) handleAnovaError(l_tmb, group, ...))
+anovaParallel <- function(list_tmb, ...) {
+  l_group <- attributes(list_tmb)$names
+  lapply(stats::setNames(l_group, l_group), function(group) handleAnovaError(list_tmb, group, ...))
 }
 
 
@@ -6441,7 +6469,7 @@ subsetGenes <- function(l_genes, mockObj) {
 ```
 
 
-```{r function-evaluationWithMixedEffect, filename =  "evaluationWithMixedEffect"}
+```{r function-evaluation_withmixedeffect, filename =  "evaluation_withmixedeffect"}
 
 #' Check if the formula contains a mixed effect structure.
 #'
@@ -6712,7 +6740,7 @@ compareInferenceToExpected <- function(tidy_tmb, ground_truth_eff, formula_used)
 
 ```
 
-```{r  test-evaluationWithMixedEffect}
+```{r  test-evaluation_withmixedeffect}
 
 
 
@@ -7205,14 +7233,14 @@ The `updateParallel()` function updates and re-fits a model for each gene. It of
 ```{r example-update, warning=FALSE,  message=FALSE}
 ## -- update your fit modifying the model family
 l_tmb <- updateParallel(formula =  kij ~ varA,
-                          l_tmb = l_tmb ,
+                          list_tmb = l_tmb ,
                           family = gaussian(), 
                           log_file = "log.txt",
                           n.cores = 1)
 
 ## -- update fit using additional model control settings
 l_tmb <- updateParallel(formula =  kij ~ varA ,
-                          l_tmb = l_tmb ,
+                          list_tmb = l_tmb ,
                           family = gaussian(), 
                           log_file = "log.txt",
                           n.cores = 1,
@@ -7222,7 +7250,7 @@ l_tmb <- updateParallel(formula =  kij ~ varA ,
 
 ## -- update your model formula and your family model
 l_tmb <- updateParallel(formula =   kij ~ varA + varB  + varA:varB ,
-                          l_tmb = l_tmb ,
+                          list_tmb = l_tmb ,
                           family = glmmTMB::nbinom2(link = "log"), 
                           log_file = "log.txt",
                           n.cores = 1)
@@ -7237,13 +7265,13 @@ Visualizing fit metrics is essential for evaluating your models. Here, we show y
 
 ```{r example-plotMetrics, warning=FALSE, message=FALSE, fig.align='center', fig.height=4, fig.width=6}
 ## -- plot all metrics
-metrics_plot(l_tmb = l_tmb)
+metrics_plot(list_tmb = l_tmb)
 ```
 
 
 ```{r example-plotMetricsFocus, warning=FALSE, message=FALSE, fig.align='center', fig.height=3, fig.width=4}
 ## -- Focus on metrics
-metrics_plot(l_tmb = l_tmb, focus = c("dispersion", "logLik"))
+metrics_plot(list_tmb = l_tmb, focus = c("dispersion", "logLik"))
 ```
 
 ## Anova to select the best model
@@ -7253,10 +7281,10 @@ Utilizing the `anovaParallel()` function enables you to perform model selection
 
 ```{r example-anova, warning=FALSE,  message=FALSE}
 ## -- update your fit modifying the model family
-l_anova <- anovaParallel(l_tmb = l_tmb)
+l_anova <- anovaParallel(list_tmb = l_tmb)
 
 ## -- additional settings
-l_anova <- anovaParallel(l_tmb = l_tmb, type = "III" )
+l_anova <- anovaParallel(list_tmb = l_tmb, type = "III" )
 
 ## -- output 
 l_anova$gene1
@@ -7416,7 +7444,7 @@ fusen::fill_description(fields = list(Title = "HTRfit"), overwrite = T)
 usethis::use_mit_license("Arnaud DUVERMY")
 usethis::use_pipe(export = TRUE)
 devtools::document()
-# Keep eval=FALSE to avoid infinite loop in case you hit the knit button
+.# 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 = "HTRfit", open_vignette = T, overwrite = T)
 ```
diff --git a/man/addBasalExpression.Rd b/man/addBasalExpression.Rd
index 12951bb8e27bac820f87e81ad13caf77de095628..2cc6230f2253efc5d3fe6bab31818470f45479f7 100644
--- a/man/addBasalExpression.Rd
+++ b/man/addBasalExpression.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/scalinggeneexpression.R
+% Please edit documentation in R/basal_expression__scaling.R
 \name{addBasalExpression}
 \alias{addBasalExpression}
 \title{Compute basal expresion for gene expression based on the coefficients data frame.}
@@ -20,7 +20,7 @@ the scaling factors for gene expression.
 \description{
 This function takes the coefficients data frame \code{dtf_coef} and computes
 basal expression for gene expression. The scaling factors are generated
-using the function \code{generate_BE}.
+using the function \code{generate_basal_expression}.
 }
 \examples{
 list_var <- init_variable()
diff --git a/man/anovaParallel.Rd b/man/anovaParallel.Rd
index 5adbe52a75f489a78a035d5db0227a5105f2b2f6..6227cbada9a5fab9100cd98d3692299094440581 100644
--- a/man/anovaParallel.Rd
+++ b/man/anovaParallel.Rd
@@ -4,10 +4,10 @@
 \alias{anovaParallel}
 \title{Perform ANOVA on Multiple glmmTMB Models in Parallel}
 \usage{
-anovaParallel(l_tmb, ...)
+anovaParallel(list_tmb, ...)
 }
 \arguments{
-\item{l_tmb}{A list of \code{glmmTMB} models, with model names corresponding to the groups.}
+\item{list_tmb}{A list of \code{glmmTMB} models, with model names corresponding to the groups.}
 
 \item{...}{Additional arguments passed to \code{\link[stats]{anova}} function.}
 }
@@ -22,7 +22,7 @@ of ANOVA results for each group.
 \examples{
 # Perform ANOVA
 data(iris)
-l_tmb<- fitModelParallel( Sepal.Length ~ Sepal.Width  + Petal.Length, 
+list_tmb<- fitModelParallel( Sepal.Length ~ Sepal.Width  + Petal.Length, 
                          data = iris, group_by = "Species", n.cores = 1 )
-anov_res <- anovaParallel(l_tmb , type = "III")
+anov_res <- anovaParallel(list_tmb , type = "III")
 }
diff --git a/man/averageByGroup.Rd b/man/averageByGroup.Rd
index 8a19b2d1fcc385c553c0836801577276a321cf8b..ec9c8ad4b0c57fb710109c9c70393725194bfb26 100644
--- a/man/averageByGroup.Rd
+++ b/man/averageByGroup.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualmainfixeffects.R
+% Please edit documentation in R/actual_mainfixeffects.R
 \name{averageByGroup}
 \alias{averageByGroup}
 \title{Calculate average values by group}
diff --git a/man/calculate_actualMixed.Rd b/man/calculate_actualMixed.Rd
index 0b907b862dafc44240d317f6f9a2c6b1ea357b89..743578921e03b4f36662aa0b41bf7aeb5a286397 100644
--- a/man/calculate_actualMixed.Rd
+++ b/man/calculate_actualMixed.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluationwithmixedeffect.R
+% Please edit documentation in R/evaluation_withmixedeffect.R
 \name{calculate_actualMixed}
 \alias{calculate_actualMixed}
 \title{Calculate actual mixed effects.}
diff --git a/man/calculate_actual_interactionX2_values.Rd b/man/calculate_actual_interactionX2_values.Rd
index 33bc332fb6ad7c1b5a1e0c8b733d2ccf8172bfa8..41ba06deb54d773680939ab4d48508e27c8066d3 100644
--- a/man/calculate_actual_interactionX2_values.Rd
+++ b/man/calculate_actual_interactionX2_values.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualinteractionfixeffects.R
+% Please edit documentation in R/actual_interactionfixeffects.R
 \name{calculate_actual_interactionX2_values}
 \alias{calculate_actual_interactionX2_values}
 \title{Calculate actual interaction values between two terms in a data frame.}
diff --git a/man/calculate_actual_interactionX3_values.Rd b/man/calculate_actual_interactionX3_values.Rd
index a5512984d07da6a694347e71b54130d0d39f9f0a..e2a7cb71adafc66527c4d1e1abf16be73f1fd129 100644
--- a/man/calculate_actual_interactionX3_values.Rd
+++ b/man/calculate_actual_interactionX3_values.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualinteractionfixeffects.R
+% Please edit documentation in R/actual_interactionfixeffects.R
 \name{calculate_actual_interactionX3_values}
 \alias{calculate_actual_interactionX3_values}
 \title{Calculate Actual Interaction Values for Three Fixed Effects}
diff --git a/man/compareInferenceToExpected.Rd b/man/compareInferenceToExpected.Rd
index 5bca8afe6a2e9783a05fa2698febabec2fd63516..da1e2bd7fe947a690c83c286bf5bbf3a048f0713 100644
--- a/man/compareInferenceToExpected.Rd
+++ b/man/compareInferenceToExpected.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluationwithmixedeffect.R
+% Please edit documentation in R/evaluation_withmixedeffect.R
 \name{compareInferenceToExpected}
 \alias{compareInferenceToExpected}
 \title{Compare inference results to expected values for a given model.}
diff --git a/man/computeActualInteractionFixEff.Rd b/man/computeActualInteractionFixEff.Rd
index 80c705fa0d9d3514aa33d4e6f9ad60f3974006b9..75dd4fcd69667590cadf26b331fb05498d6fd743 100644
--- a/man/computeActualInteractionFixEff.Rd
+++ b/man/computeActualInteractionFixEff.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualinteractionfixeffects.R
+% Please edit documentation in R/actual_interactionfixeffects.R
 \name{computeActualInteractionFixEff}
 \alias{computeActualInteractionFixEff}
 \title{Compute actual interaction values for multiple interaction terms.}
diff --git a/man/convert2Factor.Rd b/man/convert2Factor.Rd
index 73e53d12f3d6d0e05f307bdcea97e68804a73850..82b9b25b87053d1abe022e5241d7d646bf7a6744 100644
--- a/man/convert2Factor.Rd
+++ b/man/convert2Factor.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualmainfixeffects.R
+% Please edit documentation in R/utils.R
 \name{convert2Factor}
 \alias{convert2Factor}
 \title{Convert specified columns to factor}
@@ -17,3 +17,11 @@ The modified data frame with specified columns converted to factors
 \description{
 Convert specified columns to factor
 }
+\examples{
+data <- data.frame( Category1 = c("A", "B", "A", "B"),
+                     Category2 = c("X", "Y", "X", "Z"),
+                     Value = 1:4,
+                     stringsAsFactors = FALSE )
+## -- Convert columns to factors
+convert2Factor(data, columns = c("Category1", "Category2"))
+}
diff --git a/man/counts_plot.Rd b/man/counts_plot.Rd
index e95bef3c9f74099d31f0b35da2bf0b1e46c34bec..ae1c17f7161cfed66ca1eb30e5f53aa869d79261 100644
--- a/man/counts_plot.Rd
+++ b/man/counts_plot.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/countsplot.R
+% Please edit documentation in R/counts_plot.R
 \name{counts_plot}
 \alias{counts_plot}
 \title{Generate a density plot of gene counts}
diff --git a/man/dispersion_plot.Rd b/man/dispersion_plot.Rd
index df34ab348c40f702437c63d3e83cfa8805655da4..2b7af5431a202daa810e3a53f2c39cbd9eb3457f 100644
--- a/man/dispersion_plot.Rd
+++ b/man/dispersion_plot.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluatedispersion.R
+% Please edit documentation in R/evaluate_dispersion.R
 \name{dispersion_plot}
 \alias{dispersion_plot}
 \title{Dispersion Evaluation Plot}
diff --git a/man/evaluateDispersion.Rd b/man/evaluateDispersion.Rd
index 73f7cf7eac954af872b9a68d52e51743c6c20484..cfe818fbeb807f3bf0fc25bf1a2550eb18d5af65 100644
--- a/man/evaluateDispersion.Rd
+++ b/man/evaluateDispersion.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluatedispersion.R
+% Please edit documentation in R/evaluate_dispersion.R
 \name{evaluateDispersion}
 \alias{evaluateDispersion}
 \title{Evaluate Dispersion Comparison}
diff --git a/man/exportReportFile.Rd b/man/exportReportFile.Rd
index d612cde5d7c1ed69f0a0bddc5bf893db61b76d9e..83d8409568d81828297506827f6099bf9b6406b8 100644
--- a/man/exportReportFile.Rd
+++ b/man/exportReportFile.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulationreport.R
+% Please edit documentation in R/simulation_report.R
 \name{exportReportFile}
 \alias{exportReportFile}
 \title{Export the Analysis Report to a File}
diff --git a/man/extractDESeqDispersion.Rd b/man/extractDESeqDispersion.Rd
deleted file mode 100644
index b0155a8ffe83fee919f441b1c06f6ff33e122ceb..0000000000000000000000000000000000000000
--- a/man/extractDESeqDispersion.Rd
+++ /dev/null
@@ -1,23 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluatedispersion.R
-\name{extractDESeqDispersion}
-\alias{extractDESeqDispersion}
-\title{Extract DESeq2 Dispersion Values}
-\usage{
-extractDESeqDispersion(deseq_wrapped)
-}
-\arguments{
-\item{deseq_wrapped}{A DESeq2 wrapped object containing dispersion values.}
-}
-\value{
-A data frame containing inferred dispersion values.
-}
-\description{
-Extracts inferred dispersion values from a DESeq2 wrapped object.
-}
-\examples{
-\dontrun{
-dispersion_df <- extractDESeqDispersion(deseq2_object)
-print(dispersion_df)
-}
-}
diff --git a/man/extractTMBDispersion.Rd b/man/extractTMBDispersion.Rd
deleted file mode 100644
index dd586fdd3c61e90b3f36fc1368921422247f6fd2..0000000000000000000000000000000000000000
--- a/man/extractTMBDispersion.Rd
+++ /dev/null
@@ -1,23 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluatedispersion.R
-\name{extractTMBDispersion}
-\alias{extractTMBDispersion}
-\title{Extract TMB Dispersion Values}
-\usage{
-extractTMBDispersion(l_tmb)
-}
-\arguments{
-\item{l_tmb}{A TMB result object containing dispersion values.}
-}
-\value{
-A data frame containing inferred dispersion values.
-}
-\description{
-Extracts inferred dispersion values from a TMB result object.
-}
-\examples{
-\dontrun{
-dispersion_df <- extractTMBDispersion(tmb_result)
-print(dispersion_df)
-}
-}
diff --git a/man/extract_ddsDispersion.Rd b/man/extract_ddsDispersion.Rd
new file mode 100644
index 0000000000000000000000000000000000000000..e809a6cf447f9f71469afec37888ceaa688324b9
--- /dev/null
+++ b/man/extract_ddsDispersion.Rd
@@ -0,0 +1,22 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/evaluate_dispersion.R
+\name{extract_ddsDispersion}
+\alias{extract_ddsDispersion}
+\title{Extract DESeq2 Dispersion Values}
+\usage{
+extract_ddsDispersion(dds_wrapped)
+}
+\arguments{
+\item{dds_wrapped}{A DESeq2 wrapped object containing dispersion values.}
+}
+\value{
+A data frame containing inferred dispersion values.
+}
+\description{
+Extracts inferred dispersion values from a DESeq2 wrapped object.
+}
+\examples{
+\dontrun{
+dispersion_df <- extract_ddsDispersion(deseq2_object)
+}
+}
diff --git a/man/extract_tmbDispersion.Rd b/man/extract_tmbDispersion.Rd
new file mode 100644
index 0000000000000000000000000000000000000000..0e0200f5592320382cd2ededa9a3d3f42f2a6122
--- /dev/null
+++ b/man/extract_tmbDispersion.Rd
@@ -0,0 +1,22 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/evaluate_dispersion.R
+\name{extract_tmbDispersion}
+\alias{extract_tmbDispersion}
+\title{Extract TMB Dispersion Values}
+\usage{
+extract_tmbDispersion(list_tmb)
+}
+\arguments{
+\item{list_tmb}{A TMB result object containing dispersion values.}
+}
+\value{
+A data frame containing inferred dispersion values.
+}
+\description{
+Extracts inferred dispersion values from a TMB result object.
+}
+\examples{
+\dontrun{
+dispersion_df <- extract_tmbDispersion(tmb_result)
+}
+}
diff --git a/man/filter_dataframe.Rd b/man/filter_dataframe.Rd
index 5a8b3e636a530b999148760c4e2c864f3fdf80fe..ae5d68748d4192145d7d0713e82934e257989691 100644
--- a/man/filter_dataframe.Rd
+++ b/man/filter_dataframe.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualinteractionfixeffects.R
+% Please edit documentation in R/actual_interactionfixeffects.R
 \name{filter_dataframe}
 \alias{filter_dataframe}
 \title{Filter DataFrame}
diff --git a/man/findAttribute.Rd b/man/findAttribute.Rd
index 2ebfc99a6077e4d9f46fc617b502e36e7abb6027..4c2ca7313ec1c0d9e170aaef9a9884fc2dbba497 100644
--- a/man/findAttribute.Rd
+++ b/man/findAttribute.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualmainfixeffects.R
+% Please edit documentation in R/actual_mainfixeffects.R
 \name{findAttribute}
 \alias{findAttribute}
 \title{Find Attribute}
diff --git a/man/dot-fitModel.Rd b/man/fitModel.Rd
similarity index 79%
rename from man/dot-fitModel.Rd
rename to man/fitModel.Rd
index 7e392eff232e639e3e53ead62ad788db8cd9166e..a54618e34929ca42a21e868908476ea1023777e5 100644
--- a/man/dot-fitModel.Rd
+++ b/man/fitModel.Rd
@@ -1,10 +1,10 @@
 % Generated by roxygen2: do not edit by hand
 % Please edit documentation in R/fitmodel.R
-\name{.fitModel}
-\alias{.fitModel}
+\name{fitModel}
+\alias{fitModel}
 \title{Fit a model using the fitModel function.}
 \usage{
-.fitModel(formula, data, ...)
+fitModel(formula, data, ...)
 }
 \arguments{
 \item{formula}{Formula specifying the model formula}
@@ -20,5 +20,5 @@ Fitted model object or NULL if there was an error
 Fit a model using the fitModel function.
 }
 \examples{
-.fitModel(formula = mpg ~ cyl + disp, data = mtcars)
+fitModel(formula = mpg ~ cyl + disp, data = mtcars)
 }
diff --git a/man/fitUpdate.Rd b/man/fitUpdate.Rd
index 8df179fccb7748fef444c370653c54621ccb9b4a..b90494895e66357b5e49ca8170193ac11acd1b5c 100644
--- a/man/fitUpdate.Rd
+++ b/man/fitUpdate.Rd
@@ -1,13 +1,13 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/updatefitmodel.R
+% Please edit documentation in R/update_fittedmodel.R
 \name{fitUpdate}
 \alias{fitUpdate}
 \title{Fit and update a GLMNB model.}
 \usage{
-fitUpdate(glmnb_obj, formula, ...)
+fitUpdate(glmm_obj, formula, ...)
 }
 \arguments{
-\item{glmnb_obj}{A GLMNB object to be updated.}
+\item{glmm_obj}{A glmmTMB object to be updated.}
 
 \item{formula}{Formula for the updated GLMNB model.}
 
diff --git a/man/generateActualForMainFixEff.Rd b/man/generateActualForMainFixEff.Rd
index 4ccfb7f0bf1407c0b8e85be70dd9a5c895bcbf3d..f8a2b40b1f65f74914bcd15690f4e8009146d77f 100644
--- a/man/generateActualForMainFixEff.Rd
+++ b/man/generateActualForMainFixEff.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualmainfixeffects.R
+% Please edit documentation in R/actual_mainfixeffects.R
 \name{generateActualForMainFixEff}
 \alias{generateActualForMainFixEff}
 \title{Generate actual values for a given term}
diff --git a/man/generateActualInteractionX2_FixEff.Rd b/man/generateActualInteractionX2_FixEff.Rd
index 9a95303f0bc8d77a5f5761d3d8d5764059369d2b..1177c00fc242b553872367d489d84f266d70d3bd 100644
--- a/man/generateActualInteractionX2_FixEff.Rd
+++ b/man/generateActualInteractionX2_FixEff.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualinteractionfixeffects.R
+% Please edit documentation in R/actual_interactionfixeffects.R
 \name{generateActualInteractionX2_FixEff}
 \alias{generateActualInteractionX2_FixEff}
 \title{Generate actual values for the interaction fixed effect.}
diff --git a/man/generateActualInteractionX3_FixEff.Rd b/man/generateActualInteractionX3_FixEff.Rd
index d94478142c417cd3426277e6c77cbbf6b76b398a..9778e3bc9e69af494b90387e67bfbdbbbdf8074d 100644
--- a/man/generateActualInteractionX3_FixEff.Rd
+++ b/man/generateActualInteractionX3_FixEff.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualinteractionfixeffects.R
+% Please edit documentation in R/actual_interactionfixeffects.R
 \name{generateActualInteractionX3_FixEff}
 \alias{generateActualInteractionX3_FixEff}
 \title{Generate Actual Interaction Values for Three Fixed Effects}
diff --git a/man/generateCountTable.Rd b/man/generateCountTable.Rd
index 6e3182d3ae8eae08b3210ba7585a39a0d8021c00..ca61466eba77257bedaa25b27c940ba87323f662 100644
--- a/man/generateCountTable.Rd
+++ b/man/generateCountTable.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/mock-rnaseq.R
+% Please edit documentation in R/mock_rnaseq.R
 \name{generateCountTable}
 \alias{generateCountTable}
 \title{Generate count table}
diff --git a/man/generateGridCombination_fromListVar.Rd b/man/generateGridCombination_fromListVar.Rd
index b8952b720c850bacb87854a9c2e3e610e8c6f4d8..b59c6ec4db282ed7524e869f84d02ae76b4111e5 100644
--- a/man/generateGridCombination_fromListVar.Rd
+++ b/man/generateGridCombination_fromListVar.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulation_initialization.R
+% Please edit documentation in R/utils.R
 \name{generateGridCombination_fromListVar}
 \alias{generateGridCombination_fromListVar}
 \title{Get grid combination from list_var}
@@ -15,3 +15,6 @@ The grid combination between variable in list_var
 \description{
 Get grid combination from list_var
 }
+\examples{
+generateGridCombination_fromListVar(init_variable())
+}
diff --git a/man/generateReplicationMatrix.Rd b/man/generateReplicationMatrix.Rd
index b958e369116db26bc00dc917c253c4d92a7819a6..2616e88d9ab453554afdb011bf4d6b662ed8dca0 100644
--- a/man/generateReplicationMatrix.Rd
+++ b/man/generateReplicationMatrix.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/mock-rnaseq.R
+% Please edit documentation in R/mock_rnaseq.R
 \name{generateReplicationMatrix}
 \alias{generateReplicationMatrix}
 \title{Generate replication matrix}
diff --git a/man/generate_BE.Rd b/man/generate_BE.Rd
deleted file mode 100644
index 375df6eaa3b3f0286c39712b8d104f25919ad0fc..0000000000000000000000000000000000000000
--- a/man/generate_BE.Rd
+++ /dev/null
@@ -1,23 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/scalinggeneexpression.R
-\name{generate_BE}
-\alias{generate_BE}
-\title{Generate BE data.}
-\usage{
-generate_BE(n_genes, basal_expression)
-}
-\arguments{
-\item{n_genes}{The number of genes to generate BE data for.}
-
-\item{basal_expression}{a numeric vector from which sample BE for eacg genes}
-}
-\value{
-A data frame containing gene IDs, BE values
-}
-\description{
-This function generates BE data for a given number of genes, in a vector of BE values.
-}
-\examples{
-generate_BE(n_genes = 100, 10)
-
-}
diff --git a/man/generate_basal_expression.Rd b/man/generate_basal_expression.Rd
new file mode 100644
index 0000000000000000000000000000000000000000..3e68b0680d085a2df751962536c6393336594a07
--- /dev/null
+++ b/man/generate_basal_expression.Rd
@@ -0,0 +1,23 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/basal_expression__scaling.R
+\name{generate_basal_expression}
+\alias{generate_basal_expression}
+\title{Generate BE data.}
+\usage{
+generate_basal_expression(n_genes, basal_expression)
+}
+\arguments{
+\item{n_genes}{The number of genes to generate BE data for.}
+
+\item{basal_expression}{a numeric vector from which sample BE for eacg genes}
+}
+\value{
+A data frame containing gene IDs, BE values
+}
+\description{
+This function generates basal expression data for a given number of genes, in a vector of basal expression values.
+}
+\examples{
+generate_basal_expression(n_genes = 100, 10)
+
+}
diff --git a/man/getActualInteractionFixEff.Rd b/man/getActualInteractionFixEff.Rd
index 5e2652ddccaab568b676e7f56ad696108d77a34e..48ff1c4180306b2e991c1726b2e6d1d94f23520c 100644
--- a/man/getActualInteractionFixEff.Rd
+++ b/man/getActualInteractionFixEff.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualinteractionfixeffects.R
+% Please edit documentation in R/actual_interactionfixeffects.R
 \name{getActualInteractionFixEff}
 \alias{getActualInteractionFixEff}
 \title{Get the actual interaction values for a given interaction term in the data.}
diff --git a/man/getActualIntercept.Rd b/man/getActualIntercept.Rd
index 83ece5d396ec0b9804cddc4bfd167419046bc2fb..2cb37ca37a10c3680cd98626c33e16d039c5429f 100644
--- a/man/getActualIntercept.Rd
+++ b/man/getActualIntercept.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualmainfixeffects.R
+% Please edit documentation in R/actual_mainfixeffects.R
 \name{getActualIntercept}
 \alias{getActualIntercept}
 \title{Get the intercept dataframe}
diff --git a/man/getActualMainFixEff.Rd b/man/getActualMainFixEff.Rd
index 882f46d56f973faeb834251201e530fedc2f2cc4..b0409e958bcd64d31b000dff494d03547c84a04d 100644
--- a/man/getActualMainFixEff.Rd
+++ b/man/getActualMainFixEff.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualmainfixeffects.R
+% Please edit documentation in R/actual_mainfixeffects.R
 \name{getActualMainFixEff}
 \alias{getActualMainFixEff}
 \title{Get actual values for non-interaction terms}
diff --git a/man/getActualMixed_typeI.Rd b/man/getActualMixed_typeI.Rd
index 93d62f7a25859c8665c3ef8f368a9d6d9bba6ca3..e32f28d9676dd0f3a8c78a8de1bb10bd39ae8fe8 100644
--- a/man/getActualMixed_typeI.Rd
+++ b/man/getActualMixed_typeI.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluationwithmixedeffect.R
+% Please edit documentation in R/evaluation_withmixedeffect.R
 \name{getActualMixed_typeI}
 \alias{getActualMixed_typeI}
 \title{Calculate actual mixed effect values for each gene.}
diff --git a/man/getBinExpression.Rd b/man/getBinExpression.Rd
index 3aa4829d6d44225f2095d5a0ead930d808122e38..7c4ecd1566e02a6e1081140eb9258d2281a92335 100644
--- a/man/getBinExpression.Rd
+++ b/man/getBinExpression.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/scalinggeneexpression.R
+% Please edit documentation in R/basal_expression__scaling.R
 \name{getBinExpression}
 \alias{getBinExpression}
 \title{Get bin expression for a data frame.}
diff --git a/man/getCategoricalVar_inFixedEffect.Rd b/man/getCategoricalVar_inFixedEffect.Rd
index c69cc002f187d5a3110f804cd11d84df480ff0eb..ad3498dedf8de15cd8373e0091195c5e14a32623 100644
--- a/man/getCategoricalVar_inFixedEffect.Rd
+++ b/man/getCategoricalVar_inFixedEffect.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluationwithmixedeffect.R
+% Please edit documentation in R/evaluation_withmixedeffect.R
 \name{getCategoricalVar_inFixedEffect}
 \alias{getCategoricalVar_inFixedEffect}
 \title{Get the categorical variable associated with the fixed effect in a type I formula.}
diff --git a/man/dot-getColumnWithSampleID.Rd b/man/getColumnWithSampleID.Rd
similarity index 76%
rename from man/dot-getColumnWithSampleID.Rd
rename to man/getColumnWithSampleID.Rd
index bc085af6c148eebf67d00f8a478f185b6b162601..288ba8714815f8667ac4755f51762cc2064fa49a 100644
--- a/man/dot-getColumnWithSampleID.Rd
+++ b/man/getColumnWithSampleID.Rd
@@ -1,10 +1,10 @@
 % Generated by roxygen2: do not edit by hand
 % Please edit documentation in R/prepare_data2fit.R
-\name{.getColumnWithSampleID}
-\alias{.getColumnWithSampleID}
+\name{getColumnWithSampleID}
+\alias{getColumnWithSampleID}
 \title{Get column name with sampleID}
 \usage{
-.getColumnWithSampleID(dtf_countsLong, metadata)
+getColumnWithSampleID(dtf_countsLong, metadata)
 }
 \arguments{
 \item{dtf_countsLong}{Long data frame of counts}
@@ -21,5 +21,5 @@ Returns the column name in the metadata data frame that corresponds to the given
 list_var <- init_variable()
 mock_data <- mock_rnaseq(list_var, n_genes = 3, 2,2, 2)
 dtf_countLong <- countMatrix_2longDtf(mock_data$counts)
-.getColumnWithSampleID(dtf_countLong, mock_data$metadata)
+getColumnWithSampleID(dtf_countLong, mock_data$metadata)
 }
diff --git a/man/getCountsTable.Rd b/man/getCountsTable.Rd
index 035c531436076d4c6477e33ad41d2ab421e4c184..2a6537c9ad3495b85ff2cfdd3ed496b06973572f 100644
--- a/man/getCountsTable.Rd
+++ b/man/getCountsTable.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulation2.R
+% Please edit documentation in R/simulation.R
 \name{getCountsTable}
 \alias{getCountsTable}
 \title{getCountsTable}
diff --git a/man/getData2computeActualFixEffect.Rd b/man/getData2computeActualFixEffect.Rd
index 59fd6369edb753049046f94ffe4e10db62569b61..c5de18d7c77104e548ffb841e5223d7decf5f629 100644
--- a/man/getData2computeActualFixEffect.Rd
+++ b/man/getData2computeActualFixEffect.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualmainfixeffects.R
+% Please edit documentation in R/actual_mainfixeffects.R
 \name{getData2computeActualFixEffect}
 \alias{getData2computeActualFixEffect}
 \title{Get data for calculating actual values}
diff --git a/man/getDispersionComparison.Rd b/man/getDispersionComparison.Rd
index f4bd2539aad432243a629df957c4f0067af08c0d..1d763de922962a1d399f1c1b9b2d9c8f0c9727d7 100644
--- a/man/getDispersionComparison.Rd
+++ b/man/getDispersionComparison.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluatedispersion.R
+% Please edit documentation in R/evaluate_dispersion.R
 \name{getDispersionComparison}
 \alias{getDispersionComparison}
 \title{Get Dispersion Comparison}
@@ -20,6 +20,5 @@ Compares inferred dispersion values with actual dispersion values.
 \examples{
 \dontrun{
 dispersion_comparison <- getDispersionComparison(inferred_disp, actual_disp)
-print(dispersion_comparison)
 }
 }
diff --git a/man/getDispersionMatrix.Rd b/man/getDispersionMatrix.Rd
index 52e76ef990f95a4b69e1c77a5a80e638729814d3..66e03e2dba46f63aea206b1848f3cff0ea388f92 100644
--- a/man/getDispersionMatrix.Rd
+++ b/man/getDispersionMatrix.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulation2.R
+% Please edit documentation in R/simulation.R
 \name{getDispersionMatrix}
 \alias{getDispersionMatrix}
 \title{getDispersionMatrix}
diff --git a/man/getEstimate_df.Rd b/man/getEstimate_df.Rd
index 3599cb3ed2ea6abea0318b6445ccf8c5dd733c9f..398d63abd6ad6a9894c6ce037c6c3e0137b9d4c5 100644
--- a/man/getEstimate_df.Rd
+++ b/man/getEstimate_df.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/wrapperdeseq2.R
+% Please edit documentation in R/wrapper_dds.R
 \name{getEstimate_df}
 \alias{getEstimate_df}
 \title{Extract Inferred Estimate Information from DESeq2 Results}
diff --git a/man/getGivenAttribute.Rd b/man/getGivenAttribute.Rd
index f29d9e522b73dac49bc0eae924373291686d7b83..f765cb9bfac76e73381e5a887cbbe1fbceb1db96 100644
--- a/man/getGivenAttribute.Rd
+++ b/man/getGivenAttribute.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulation_initialization.R
+% Please edit documentation in R/utils.R
 \name{getGivenAttribute}
 \alias{getGivenAttribute}
 \title{Get a given attribute from a list of variables}
@@ -7,7 +7,7 @@
 getGivenAttribute(list_var, attribute)
 }
 \arguments{
-\item{list_var}{A list of variables (already initialized)}
+\item{list_var}{A list of variables (already initialized with init_variable)}
 
 \item{attribute}{A string specifying the attribute to retrieve in all occurrences of the list}
 }
@@ -17,3 +17,6 @@ A list without NULL values
 \description{
 Get a given attribute from a list of variables
 }
+\examples{
+getGivenAttribute(init_variable(), "level")
+}
diff --git a/man/getGlance.Rd b/man/getGlance.Rd
index 4b79867df3ced684dc8b4c57c1118c8f3f1fd1e1..239baf8810cee85108b8147d5c6804f4f06f9de2 100644
--- a/man/getGlance.Rd
+++ b/man/getGlance.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/glance_tmb.R
+% Please edit documentation in R/glance_glmmtmb.R
 \name{getGlance}
 \alias{getGlance}
 \title{Extracts the summary statistics from a single glmmTMB model.}
diff --git a/man/getGridCombination.Rd b/man/getGridCombination.Rd
index b5489f6282f89587597203d8b62ead00f3bd97d8..4617d767183a190d27eda1a6f4c8903e0f26ed98 100644
--- a/man/getGridCombination.Rd
+++ b/man/getGridCombination.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulation_initialization.R
+% Please edit documentation in R/utils.R
 \name{getGridCombination}
 \alias{getGridCombination}
 \title{getGridCombination}
diff --git a/man/getGrobTable.Rd b/man/getGrobTable.Rd
index f15c04d5e4edacea39a1c4c05b930ba3e71781bd..2fe06cb85b6f18774b827d9e320b7cae0123538b 100644
--- a/man/getGrobTable.Rd
+++ b/man/getGrobTable.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulationreport.R
+% Please edit documentation in R/simulation_report.R
 \name{getGrobTable}
 \alias{getGrobTable}
 \title{Generate a Formatted Table as a Grid Graphics Object}
diff --git a/man/getLabelExpected.Rd b/man/getLabelExpected.Rd
index 0b9791f0490dbf02e834c697e947003bd4a21801..0a1147e4a05a55c325946eed9de9ae844d4626bd 100644
--- a/man/getLabelExpected.Rd
+++ b/man/getLabelExpected.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/rocplot.R
+% Please edit documentation in R/roc_plot.R
 \name{getLabelExpected}
 \alias{getLabelExpected}
 \title{Get Labels for Expected Differential Expression}
diff --git a/man/getLabels.Rd b/man/getLabels.Rd
index df259d30d33b010b23d2090a9a2f77d0de3725fb..a89dd51d26c5543c4acacb798204c173b60321ce 100644
--- a/man/getLabels.Rd
+++ b/man/getLabels.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulation_initialization.R
+% Please edit documentation in R/utils.R
 \name{getLabels}
 \alias{getLabels}
 \title{Get labels for variables}
@@ -17,3 +17,6 @@ A list of labels per variable
 \description{
 Get labels for variables
 }
+\examples{
+labels <- getLabels(c("varA", "varB"), c(2, 3))
+}
diff --git a/man/getListVar.Rd b/man/getListVar.Rd
index 88257b2624249ec94b77119fa56b0c047825a766..b67c882b9f97402464739db8374fe644ebc3b9a6 100644
--- a/man/getListVar.Rd
+++ b/man/getListVar.Rd
@@ -1,13 +1,13 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulation_initialization.R
+% Please edit documentation in R/utils.R
 \name{getListVar}
 \alias{getListVar}
 \title{Get the list of variable names}
 \usage{
-getListVar(input)
+getListVar(list_var)
 }
 \arguments{
-\item{input}{R list, e.g., output of init_variable}
+\item{list_var}{R list, e.g., output of init_variable}
 }
 \value{
 A character vector with the names of variables
@@ -15,3 +15,6 @@ A character vector with the names of variables
 \description{
 Get the list of variable names
 }
+\examples{
+getListVar(init_variable())
+}
diff --git a/man/getLog_qij.Rd b/man/getLog_qij.Rd
index e5ce30cacbb3dc324c2156660bac0a5712897ccf..8e580d5791d9ac222911fba00ba7c079386ce6a1 100644
--- a/man/getLog_qij.Rd
+++ b/man/getLog_qij.Rd
@@ -15,3 +15,8 @@ The coefficient data frame with log_qij column added.
 \description{
 Get the log_qij values from the coefficient data frame.
 }
+\examples{
+list_var <- init_variable()
+dtf_coef <- getInput2simulation(list_var, 10)
+dtf_coef <- getLog_qij(dtf_coef)
+}
diff --git a/man/getMu_ij_matrix.Rd b/man/getMu_ij_matrix.Rd
index 67a057619ba159042f4b7b836688348a46299e75..a8f31cf99369287ecc12d6b6e396a75c846f8058 100644
--- a/man/getMu_ij_matrix.Rd
+++ b/man/getMu_ij_matrix.Rd
@@ -15,3 +15,11 @@ A Mu_ij matrix.
 \description{
 Get the Mu_ij matrix.
 }
+\examples{
+list_var <- init_variable()
+dtf_coef <- getInput2simulation(list_var, 10)
+dtf_coef <- getLog_qij(dtf_coef)
+dtf_coef <- addBasalExpression(dtf_coef, 10, c(10, 20, 0))
+dtf_coef<- getMu_ij(dtf_coef)
+getMu_ij_matrix(dtf_coef)
+}
diff --git a/man/getReplicationMatrix.Rd b/man/getReplicationMatrix.Rd
index e0edd1bdb39f21909ce477497b4439bf0bfa9ef5..84aa9847077b8c8d67893835dd9ede266909160b 100644
--- a/man/getReplicationMatrix.Rd
+++ b/man/getReplicationMatrix.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulation2.R
+% Please edit documentation in R/simulation.R
 \name{getReplicationMatrix}
 \alias{getReplicationMatrix}
 \title{getReplicationMatrix}
diff --git a/man/getSE_df.Rd b/man/getSE_df.Rd
index 61533dddcece67821284bf599176bace260975f7..37be566a2583477726fbf0e44f250610c2fc03a7 100644
--- a/man/getSE_df.Rd
+++ b/man/getSE_df.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/wrapperdeseq2.R
+% Please edit documentation in R/wrapper_dds.R
 \name{getSE_df}
 \alias{getSE_df}
 \title{Extract Standard Error Information from DESeq2 Results}
diff --git a/man/getSampleID.Rd b/man/getSampleID.Rd
index 81d4301faaf93291e9da17458ee73e1e73c824ff..4f9ab1daa29680e3b9d8a415710138e21c014547 100644
--- a/man/getSampleID.Rd
+++ b/man/getSampleID.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulation2.R
+% Please edit documentation in R/simulation.R
 \name{getSampleID}
 \alias{getSampleID}
 \title{getSampleID}
@@ -15,3 +15,6 @@ A sorted vector of sample IDs
 \description{
 getSampleID
 }
+\examples{
+getSampleID(init_variable())
+}
diff --git a/man/getSampleMetadata.Rd b/man/getSampleMetadata.Rd
index f16657224d2c2cf242b2059d5e1ef0a8d255f97e..b32827af1d43a03c9813cf2dc556cba07621e38c 100644
--- a/man/getSampleMetadata.Rd
+++ b/man/getSampleMetadata.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulation2.R
+% Please edit documentation in R/simulation.R
 \name{getSampleMetadata}
 \alias{getSampleMetadata}
 \title{Get sample metadata}
diff --git a/man/getValidDispersion.Rd b/man/getValidDispersion.Rd
index aac2e249bfd1246654b42ef4ad272d483a4289cc..194847eacf9316e8190b7ba03db1c18ad2cc1d43 100644
--- a/man/getValidDispersion.Rd
+++ b/man/getValidDispersion.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/mock-rnaseq.R
+% Please edit documentation in R/mock_rnaseq.R
 \name{getValidDispersion}
 \alias{getValidDispersion}
 \title{Validate and Filter Dispersion Values}
diff --git a/man/get_inference.Rd b/man/get_inference_dds.Rd
similarity index 79%
rename from man/get_inference.Rd
rename to man/get_inference_dds.Rd
index 3d7864a6dcd7e3e07a3ba766942ab8e42e58a3cc..e8e9cb16d27109859d7720121a3f9e0e514d2e5f 100644
--- a/man/get_inference.Rd
+++ b/man/get_inference_dds.Rd
@@ -1,10 +1,10 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/wrapperdeseq2.R
-\name{get_inference}
-\alias{get_inference}
+% Please edit documentation in R/wrapper_dds.R
+\name{get_inference_dds}
+\alias{get_inference_dds}
 \title{Calculate Inference for Differential Expression Analysis}
 \usage{
-get_inference(dds_full, lfcThreshold, altHypothesis, correction_method)
+get_inference_dds(dds_full, lfcThreshold, altHypothesis, correction_method)
 }
 \arguments{
 \item{dds_full}{A data frame containing DESeq2 results, including estimate and standard error information.}
@@ -24,7 +24,7 @@ This function calculates inference for differential expression analysis based on
 \examples{
 \dontrun{
 # Example usage of the function
-inference_result <- get_inference(dds_full, lfcThreshold = 0.5, 
+inference_result <- get_inference_dds(dds_full, lfcThreshold = 0.5, 
                                    altHypothesis = "greater", 
                                    correction_method = "BH")
 }
diff --git a/man/glance_tmb.Rd b/man/glance_tmb.Rd
index 9f57dea776ec460ee8b145293f33fe0c258b9b32..276f7783cc7dc5c5e043951b9b93dd689d158f27 100644
--- a/man/glance_tmb.Rd
+++ b/man/glance_tmb.Rd
@@ -1,13 +1,13 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/glance_tmb.R
+% Please edit documentation in R/glance_glmmtmb.R
 \name{glance_tmb}
 \alias{glance_tmb}
 \title{Extracts the summary statistics from a list of glmmTMB models.}
 \usage{
-glance_tmb(l_tmb)
+glance_tmb(list_tmb)
 }
 \arguments{
-\item{l_tmb}{A list of glmmTMB models or a unique glmmTMB obj model}
+\item{list_tmb}{A list of glmmTMB models or a unique glmmTMB obj model}
 }
 \value{
 A DataFrame with the summary statistics for all the glmmTMB models in the list.
diff --git a/man/group_logQij_per_genes_and_labels.Rd b/man/group_logQij_per_genes_and_labels.Rd
index 05aa7d7e21e9827e66b74f241a18f5068912c789..ed3f3176c03a6b05e741caf76f7f6be9005274c2 100644
--- a/man/group_logQij_per_genes_and_labels.Rd
+++ b/man/group_logQij_per_genes_and_labels.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluationwithmixedeffect.R
+% Please edit documentation in R/evaluation_withmixedeffect.R
 \name{group_logQij_per_genes_and_labels}
 \alias{group_logQij_per_genes_and_labels}
 \title{Group log_qij values per genes and labels.}
diff --git a/man/handleAnovaError.Rd b/man/handleAnovaError.Rd
index 156749c796cf6dcaf247c56ecc260a70d96fdcff..da2daa6db97161cab87bf1d180ce70dd76cda39c 100644
--- a/man/handleAnovaError.Rd
+++ b/man/handleAnovaError.Rd
@@ -4,10 +4,10 @@
 \alias{handleAnovaError}
 \title{Handle ANOVA Errors}
 \usage{
-handleAnovaError(l_TMB, group, ...)
+handleAnovaError(list_tmb, group, ...)
 }
 \arguments{
-\item{l_TMB}{A list of fitted glmmTMB models.}
+\item{list_tmb}{A list of fitted glmmTMB models.}
 
 \item{group}{A character string indicating the group for which ANOVA is calculated.}
 
@@ -20,8 +20,8 @@ A data frame containing ANOVA results for the specified group.
 This function handles ANOVA errors and warnings during the ANOVA calculation process.
 }
 \examples{
-l_tmb <- fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length,
+list_tmb <- fitModelParallel(Sepal.Length ~ Sepal.Width + Petal.Length,
                           data = iris, group_by = "Species", n.cores = 1)
-anova_res <- handleAnovaError(l_tmb, "setosa", type = "III")
+anova_res <- handleAnovaError(list_tmb, "setosa", type = "III")
 
 }
diff --git a/man/identity_plot.Rd b/man/identity_plot.Rd
index 8d1223e1fb754b7ddd649a42b2b340932d564fe2..858b776afac9e64884f4d58a98ff37d111640144 100644
--- a/man/identity_plot.Rd
+++ b/man/identity_plot.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/identityplot.R
+% Please edit documentation in R/identity_plot.R
 \name{identity_plot}
 \alias{identity_plot}
 \title{Generate an identity plot}
diff --git a/man/inferenceToExpected_withMixedEff.Rd b/man/inferenceToExpected_withMixedEff.Rd
index 5e177211789ff923eb1bec9465cbe01e127db519..bcb6dff787af2e14d8408b95cad8bf417992be49 100644
--- a/man/inferenceToExpected_withMixedEff.Rd
+++ b/man/inferenceToExpected_withMixedEff.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluationwithmixedeffect.R
+% Please edit documentation in R/evaluation_withmixedeffect.R
 \name{inferenceToExpected_withMixedEff}
 \alias{inferenceToExpected_withMixedEff}
 \title{Compare the mixed-effects inference to expected values.}
diff --git a/man/isValidInput2fit.Rd b/man/isValidInput2fit.Rd
index 40eb9278af6cb0e795d3dfe9dcf62ed5208d6af8..368d43d5d9c081af485df89a90800a08ad543773 100644
--- a/man/isValidInput2fit.Rd
+++ b/man/isValidInput2fit.Rd
@@ -22,4 +22,3 @@ data(iris)
 formula <- Sepal.Length ~ Sepal.Width + Petal.Length
 isValidInput2fit(iris, formula) # Returns TRUE if all required variables are present
 }
-\keyword{internal}
diff --git a/man/dot-isDispersionMatrixValid.Rd b/man/is_dispersionMatrixValid.Rd
similarity index 66%
rename from man/dot-isDispersionMatrixValid.Rd
rename to man/is_dispersionMatrixValid.Rd
index dcaaae9e4198eaee49a134e51e3092c1b2fa0cf9..e84a486acb76a8a67100e96f7857591e15d19c60 100644
--- a/man/dot-isDispersionMatrixValid.Rd
+++ b/man/is_dispersionMatrixValid.Rd
@@ -1,10 +1,10 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/mock-rnaseq.R
-\name{.isDispersionMatrixValid}
-\alias{.isDispersionMatrixValid}
+% Please edit documentation in R/mock_rnaseq.R
+\name{is_dispersionMatrixValid}
+\alias{is_dispersionMatrixValid}
 \title{Check the validity of the dispersion matrix}
 \usage{
-.isDispersionMatrixValid(matx_dispersion, matx_bool_replication)
+is_dispersionMatrixValid(matx_dispersion, matx_bool_replication)
 }
 \arguments{
 \item{matx_dispersion}{Replication matrix}
@@ -20,5 +20,5 @@ Checks if the dispersion matrix has the correct dimensions.
 \examples{
 matx_dispersion <- matrix(1:12, nrow = 3, ncol = 4)
 matx_bool_replication <- matrix(TRUE, nrow = 3, ncol = 4)
-.isDispersionMatrixValid(matx_dispersion, matx_bool_replication)
+is_dispersionMatrixValid(matx_dispersion, matx_bool_replication)
 }
diff --git a/man/is_formula_mixedTypeI.Rd b/man/is_formula_mixedTypeI.Rd
index 9ee9c985f2e19324b08a155d94da901863089118..80c8949392f3a9685fbf2dd97cd537f7f2a58395 100644
--- a/man/is_formula_mixedTypeI.Rd
+++ b/man/is_formula_mixedTypeI.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluationwithmixedeffect.R
+% Please edit documentation in R/evaluation_withmixedeffect.R
 \name{is_formula_mixedTypeI}
 \alias{is_formula_mixedTypeI}
 \title{Check if the formula follows a specific type I mixed effect structure.}
diff --git a/man/is_mixedEffect_inFormula.Rd b/man/is_mixedEffect_inFormula.Rd
index 69f0e17117229518f0ed772fa7c160f78066d082..8e6eca6b1fb0eeab0b5bbf2875f226249efd9b7e 100644
--- a/man/is_mixedEffect_inFormula.Rd
+++ b/man/is_mixedEffect_inFormula.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/evaluationwithmixedeffect.R
+% Please edit documentation in R/evaluation_withmixedeffect.R
 \name{is_mixedEffect_inFormula}
 \alias{is_mixedEffect_inFormula}
 \title{Check if the formula contains a mixed effect structure.}
diff --git a/man/is_positive_definite.Rd b/man/is_positive_definite.Rd
index 46783066ef744cd900d22991ff7bf01f299b10de..b96c1db488724f422ce09d2a05be68ee05031e02 100644
--- a/man/is_positive_definite.Rd
+++ b/man/is_positive_definite.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/datafrommvrnorm_manipulations.R
+% Please edit documentation in R/utils.R
 \name{is_positive_definite}
 \alias{is_positive_definite}
 \title{Check if a matrix is positive definite
diff --git a/man/launchUpdate.Rd b/man/launchUpdate.Rd
index cb9a8796b5ffd93df1b7f4887946f4b6e013d89e..7aa0c4c234e9df760c702f56e57309062d16a9ee 100644
--- a/man/launchUpdate.Rd
+++ b/man/launchUpdate.Rd
@@ -1,13 +1,13 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/updatefitmodel.R
+% Please edit documentation in R/update_fittedmodel.R
 \name{launchUpdate}
 \alias{launchUpdate}
 \title{Launch the update process for a GLMNB model.}
 \usage{
-launchUpdate(glmnb_obj, formula, ...)
+launchUpdate(glmm_obj, formula, ...)
 }
 \arguments{
-\item{glmnb_obj}{A GLMNB object to be updated.}
+\item{glmm_obj}{A glmmTMB object to be updated.}
 
 \item{formula}{Formula for the updated GLMNB model.}
 
diff --git a/man/metrics_plot.Rd b/man/metrics_plot.Rd
index f01ebf10c22d090eb53bf5e473ecbdc2af45437f..89a21126074d2f9c99cead64adf82478cf2d1008 100644
--- a/man/metrics_plot.Rd
+++ b/man/metrics_plot.Rd
@@ -4,10 +4,10 @@
 \alias{metrics_plot}
 \title{Plot Metrics for Generalized Linear Mixed Models (GLMM)}
 \usage{
-metrics_plot(l_tmb, focus = NULL)
+metrics_plot(list_tmb, focus = NULL)
 }
 \arguments{
-\item{l_tmb}{A list of GLMM objects to extract metrics from.}
+\item{list_tmb}{A list of GLMM objects to extract metrics from.}
 
 \item{focus}{A character vector specifying the metrics to focus on. Possible
 values include "AIC", "BIC", "logLik", "deviance", "df.resid", and
diff --git a/man/mock_rnaseq.Rd b/man/mock_rnaseq.Rd
index 6f67997bdeba80905ec1e111990424e0b068a3a6..e960338a6aed97cf713451def9271137f22b8f04 100644
--- a/man/mock_rnaseq.Rd
+++ b/man/mock_rnaseq.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/mock-rnaseq.R
+% Please edit documentation in R/mock_rnaseq.R
 \name{mock_rnaseq}
 \alias{mock_rnaseq}
 \title{Perform RNA-seq simulation}
diff --git a/man/dot-parallel_fit.Rd b/man/parallel_fit.Rd
similarity index 86%
rename from man/dot-parallel_fit.Rd
rename to man/parallel_fit.Rd
index af8d9d512366fa0e6f9d24e3683f8d1c774fdf9e..0474bfb76db7c5270f829ce127b8e86d2e90e9e1 100644
--- a/man/dot-parallel_fit.Rd
+++ b/man/parallel_fit.Rd
@@ -1,11 +1,11 @@
 % Generated by roxygen2: do not edit by hand
 % Please edit documentation in R/fitmodel.R
-\name{.parallel_fit}
-\alias{.parallel_fit}
+\name{parallel_fit}
+\alias{parallel_fit}
 \title{Fit models in parallel for each group using mclapply and handle logging.
 Uses parallel_fit to fit the models.}
 \usage{
-.parallel_fit(groups, group_by, formula, data, n.cores = NULL, log_file, ...)
+parallel_fit(groups, group_by, formula, data, n.cores = NULL, log_file, ...)
 }
 \arguments{
 \item{groups}{Vector of unique group values}
@@ -31,7 +31,7 @@ Fit models in parallel for each group using mclapply and handle logging.
 Uses parallel_fit to fit the models.
 }
 \examples{
-.parallel_fit(group_by = "Species", "setosa", 
+parallel_fit(group_by = "Species", "setosa", 
                formula = Sepal.Length ~ Sepal.Width + Petal.Length, 
                data = iris, n.cores = 1, log_file = "log.txt" )
 }
diff --git a/man/dot-parallel_update.Rd b/man/parallel_update.Rd
similarity index 63%
rename from man/dot-parallel_update.Rd
rename to man/parallel_update.Rd
index f3a251b3992f048966ae785f9d917887b95f08a9..ca1c208b604618dbc04b45e9917598643e83d4b2 100644
--- a/man/dot-parallel_update.Rd
+++ b/man/parallel_update.Rd
@@ -1,15 +1,15 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/updatefitmodel.R
-\name{.parallel_update}
-\alias{.parallel_update}
-\title{Internal function to fit GLMNB models in parallel.}
+% Please edit documentation in R/update_fittedmodel.R
+\name{parallel_update}
+\alias{parallel_update}
+\title{Internal function to fit glmmTMB models in parallel.}
 \usage{
-.parallel_update(formula, l_tmb, n.cores = NULL, log_file = "log.txt", ...)
+parallel_update(formula, list_tmb, n.cores = NULL, log_file = "log.txt", ...)
 }
 \arguments{
 \item{formula}{Formula for the GLMNB model.}
 
-\item{l_tmb}{List of GLMNB objects.}
+\item{list_tmb}{List of glmmTMB objects.}
 
 \item{n.cores}{Number of cores to use for parallel processing.}
 
@@ -21,7 +21,7 @@
 A list of updated GLMNB models.
 }
 \description{
-This function is used internally by \code{\link{updateParallel}} to fit GLMNB models in parallel.
+This function is used internally by \code{\link{updateParallel}} to fit glmmTMB models in parallel.
 }
 \examples{
 data(iris)
@@ -30,5 +30,5 @@ group_by <- "Species"
 formula <- Sepal.Length ~ Sepal.Width + Petal.Length
 fitted_models <- fitModelParallel(formula, iris, group_by, n.cores = 1)
 new_formula <- Sepal.Length ~ Sepal.Width 
-results <- .parallel_update(new_formula, fitted_models, n.cores = 1)
+results <- parallel_update(new_formula, fitted_models, n.cores = 1)
 }
diff --git a/man/prepareData2computeInteraction.Rd b/man/prepareData2computeInteraction.Rd
index c7e81c729ea64748658ad74eaa4fbdfc7aa9b989..9fce5ca5577db0aaaaf66e342db5b6e5b8c6164a 100644
--- a/man/prepareData2computeInteraction.Rd
+++ b/man/prepareData2computeInteraction.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualinteractionfixeffects.R
+% Please edit documentation in R/actual_interactionfixeffects.R
 \name{prepareData2computeInteraction}
 \alias{prepareData2computeInteraction}
 \title{Prepare data for computing interaction values.}
diff --git a/man/removeDuplicatedWord.Rd b/man/removeDuplicatedWord.Rd
index 9b0938b56417f2d90000bdb0cd17f76bdcb5626c..1e1e780101bdcf8a4be4a70dae76cc6f2d497c22 100644
--- a/man/removeDuplicatedWord.Rd
+++ b/man/removeDuplicatedWord.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/tidy_glmmtmb.R
+% Please edit documentation in R/utils.R
 \name{removeDuplicatedWord}
 \alias{removeDuplicatedWord}
 \title{Remove Duplicated Words from Strings}
@@ -16,7 +16,6 @@ A character vector with duplicated words removed from each string.
 This function takes a vector of strings and removes duplicated words within each string.
 }
 \examples{
-
 words <- c("hellohello", "worldworld", "programmingprogramming", "R isis great")
 cleaned_words <- removeDuplicatedWord(words)
 }
diff --git a/man/reorderColumns.Rd b/man/reorderColumns.Rd
index 1204e33636cc803f4a38ec0b3c58599bd74e927e..0e1cf11b1e808d8aa5d6b85e822a0679fc89b158 100644
--- a/man/reorderColumns.Rd
+++ b/man/reorderColumns.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/tidy_glmmtmb.R
+% Please edit documentation in R/utils.R
 \name{reorderColumns}
 \alias{reorderColumns}
 \title{Reorder the columns of a dataframe}
diff --git a/man/dot-replicateByGroup.Rd b/man/replicateByGroup.Rd
similarity index 73%
rename from man/dot-replicateByGroup.Rd
rename to man/replicateByGroup.Rd
index c50c45a84f6b7c1518675edfd7fa90d433ce8e94..719d7bef1eabca415eb9995d710fdc6532014f5e 100644
--- a/man/dot-replicateByGroup.Rd
+++ b/man/replicateByGroup.Rd
@@ -1,10 +1,10 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulation2.R
-\name{.replicateByGroup}
-\alias{.replicateByGroup}
+% Please edit documentation in R/simulation.R
+\name{replicateByGroup}
+\alias{replicateByGroup}
 \title{Replicate rows of a data frame by group}
 \usage{
-.replicateByGroup(df, group_var, rep_list)
+replicateByGroup(df, group_var, rep_list)
 }
 \arguments{
 \item{df}{Data frame to replicate}
@@ -21,6 +21,6 @@ Replicates the rows of a data frame based on a grouping variable and replication
 }
 \examples{
 df <- data.frame(group = c("A", "B"), value = c(1, 2))
-.replicateByGroup(df, "group", c(2, 3))
+replicateByGroup(df, "group", c(2, 3))
 
 }
diff --git a/man/dot-replicateMatrix.Rd b/man/replicateMatrix.Rd
similarity index 67%
rename from man/dot-replicateMatrix.Rd
rename to man/replicateMatrix.Rd
index 3e775adc84d219f9e7ecefbdf0a3ce9de52686f8..0e276e0ca01b51a531401d62e7969af5676bbf8e 100644
--- a/man/dot-replicateMatrix.Rd
+++ b/man/replicateMatrix.Rd
@@ -1,10 +1,10 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/mock-rnaseq.R
-\name{.replicateMatrix}
-\alias{.replicateMatrix}
+% Please edit documentation in R/mock_rnaseq.R
+\name{replicateMatrix}
+\alias{replicateMatrix}
 \title{Replicate matrix}
 \usage{
-.replicateMatrix(matrix, replication_matrix)
+replicateMatrix(matrix, replication_matrix)
 }
 \arguments{
 \item{matrix}{Matrix to replicate}
@@ -20,5 +20,5 @@ Replicates a matrix based on a replication matrix.
 \examples{
 matrix <- matrix(1:9, nrow = 3, ncol = 3)
 replication_matrix <- matrix(TRUE, nrow = 3, ncol = 3)
-.replicateMatrix(matrix, replication_matrix)
+replicateMatrix(matrix, replication_matrix)
 }
diff --git a/man/dot-replicateRows.Rd b/man/replicateRows.Rd
similarity index 72%
rename from man/dot-replicateRows.Rd
rename to man/replicateRows.Rd
index 497910b80a40291e9c1f73be1d2f5ed9d4592e46..4eef78001d9ba005e426fb45d69f1d8c176749cf 100644
--- a/man/dot-replicateRows.Rd
+++ b/man/replicateRows.Rd
@@ -1,10 +1,10 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulation2.R
-\name{.replicateRows}
-\alias{.replicateRows}
+% Please edit documentation in R/simulation.R
+\name{replicateRows}
+\alias{replicateRows}
 \title{Replicate rows of a data frame}
 \usage{
-.replicateRows(df, n)
+replicateRows(df, n)
 }
 \arguments{
 \item{df}{Data frame to replicate}
@@ -19,6 +19,6 @@ Replicates the rows of a data frame by a specified factor.
 }
 \examples{
 df <- data.frame(a = 1:3, b = letters[1:3])
-.replicateRows(df, 2)
+replicateRows(df, 2)
 
 }
diff --git a/man/roc_plot.Rd b/man/roc_plot.Rd
index 7b6ed205c65fe04b461bf65354374c16a1b44dcb..bb68dfef407cfdcbb1cfb6705dc288ed69d4ce80 100644
--- a/man/roc_plot.Rd
+++ b/man/roc_plot.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/rocplot.R
+% Please edit documentation in R/roc_plot.R
 \name{roc_plot}
 \alias{roc_plot}
 \title{Generate ROC Curve Plot}
diff --git a/man/scaleCountsTable.Rd b/man/scaleCountsTable.Rd
index 349f87cc537d585271551aeef32d3001ac47ba52..1ad7280869fd056437cdfa9c1b4bff4332a30964 100644
--- a/man/scaleCountsTable.Rd
+++ b/man/scaleCountsTable.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/scalingsequencingdepth.R
+% Please edit documentation in R/sequencing_depth_scaling.R
 \name{scaleCountsTable}
 \alias{scaleCountsTable}
 \title{Scale Counts Table}
diff --git a/man/simulationReport.Rd b/man/simulationReport.Rd
index 123a3af603ad14486e6d3fd7bf943dceb68090cc..9f4a148867c2145ba38103826fbe8218c3ba0873 100644
--- a/man/simulationReport.Rd
+++ b/man/simulationReport.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/simulationreport.R
+% Please edit documentation in R/simulation_report.R
 \name{simulationReport}
 \alias{simulationReport}
 \title{Generate a simulation report}
diff --git a/man/subsetByTermLabel.Rd b/man/subsetByTermLabel.Rd
index 95138fb8f4a4ab4e99dde1f7418177c1c3a333d0..9a0d0ff175e85fb99e515bdf01ae3304094bd76b 100644
--- a/man/subsetByTermLabel.Rd
+++ b/man/subsetByTermLabel.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualmainfixeffects.R
+% Please edit documentation in R/actual_mainfixeffects.R
 \name{subsetByTermLabel}
 \alias{subsetByTermLabel}
 \title{subset data By Term Label}
diff --git a/man/dot-subsetData_andfit.Rd b/man/subsetData_andfit.Rd
similarity index 80%
rename from man/dot-subsetData_andfit.Rd
rename to man/subsetData_andfit.Rd
index f7db544f479ad270c0a668e1c1efd42519afe7ee..e84b1baf65417fc4723a2ebb6f294ebcb6dc1529 100644
--- a/man/dot-subsetData_andfit.Rd
+++ b/man/subsetData_andfit.Rd
@@ -1,10 +1,10 @@
 % Generated by roxygen2: do not edit by hand
 % Please edit documentation in R/fitmodel.R
-\name{.subsetData_andfit}
-\alias{.subsetData_andfit}
+\name{subsetData_andfit}
+\alias{subsetData_andfit}
 \title{Fit the model based using fitModel functions.}
 \usage{
-.subsetData_andfit(group, group_by, formula, data, ...)
+subsetData_andfit(group, group_by, formula, data, ...)
 }
 \arguments{
 \item{group}{The specific group to fit the model for}
@@ -24,7 +24,7 @@ Fitted model object or NULL if there was an error
 Fit the model based using fitModel functions.
 }
 \examples{
-.subsetData_andfit(group = "setosa", group_by = "Species", 
+subsetData_andfit(group = "setosa", group_by = "Species", 
                  formula = Sepal.Length ~ Sepal.Width + Petal.Length, 
                  data = iris )
 }
diff --git a/man/subsetFixEffectInferred.Rd b/man/subsetFixEffectInferred.Rd
index aebd49bec3e14408f0b7317d35b1a84942701690..a6c96ec8574e282498c6aff84f6cd0fb7d16324d 100644
--- a/man/subsetFixEffectInferred.Rd
+++ b/man/subsetFixEffectInferred.Rd
@@ -1,5 +1,5 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/actualmainfixeffects.R
+% Please edit documentation in R/actual_mainfixeffects.R
 \name{subsetFixEffectInferred}
 \alias{subsetFixEffectInferred}
 \title{Subset Fixed Effect Inferred Terms}
diff --git a/man/tidy_tmb.Rd b/man/tidy_tmb.Rd
index 34688772aadd23ec9dc6dbcb95f456aa44832d2f..3e441ad3563c71b881d132b0f98126030b27857a 100644
--- a/man/tidy_tmb.Rd
+++ b/man/tidy_tmb.Rd
@@ -4,10 +4,10 @@
 \alias{tidy_tmb}
 \title{Extract Tidy Summary of Multiple glmmTMB Models}
 \usage{
-tidy_tmb(l_tmb)
+tidy_tmb(list_tmb)
 }
 \arguments{
-\item{l_tmb}{A list of glmmTMB model objects.}
+\item{list_tmb}{A list of glmmTMB model objects.}
 }
 \value{
 A data frame containing a tidy summary of the fixed and random effects from all glmmTMB models in the list.
diff --git a/man/updateParallel.Rd b/man/updateParallel.Rd
index ca48ac4da540fcf572b734efa792e8b5a1d2d410..4c0b838a1199cb7045754dc7cac33e6e5e84e633 100644
--- a/man/updateParallel.Rd
+++ b/man/updateParallel.Rd
@@ -1,15 +1,15 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/updatefitmodel.R
+% Please edit documentation in R/update_fittedmodel.R
 \name{updateParallel}
 \alias{updateParallel}
-\title{Update GLMNB models in parallel.}
+\title{Update glmmTMB models in parallel.}
 \usage{
-updateParallel(formula, l_tmb, n.cores = NULL, log_file = "log.txt", ...)
+updateParallel(formula, list_tmb, n.cores = NULL, log_file = "log.txt", ...)
 }
 \arguments{
 \item{formula}{Formula for the GLMNB model.}
 
-\item{l_tmb}{List of GLMNB objects.}
+\item{list_tmb}{List of glmmTMB objects.}
 
 \item{n.cores}{Number of cores to use for parallel processing. If NULL, the function will use all available cores.}
 
@@ -21,7 +21,7 @@ updateParallel(formula, l_tmb, n.cores = NULL, log_file = "log.txt", ...)
 A list of updated GLMNB models.
 }
 \description{
-This function fits GLMNB models in parallel using multiple cores, allowing for faster computation.
+This function fits glmmTMB models in parallel using multiple cores, allowing for faster computation.
 }
 \examples{
 data(iris)
diff --git a/man/wrapper_DESeq2.Rd b/man/wrap_dds.Rd
similarity index 87%
rename from man/wrapper_DESeq2.Rd
rename to man/wrap_dds.Rd
index 2b3f5b388c8316e85b1f8568cab43dcbe026c164..5d798d4cb57a5a42d4ca17e946d735a701cababd 100644
--- a/man/wrapper_DESeq2.Rd
+++ b/man/wrap_dds.Rd
@@ -1,10 +1,10 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/wrapperdeseq2.R
-\name{wrapper_DESeq2}
-\alias{wrapper_DESeq2}
+% Please edit documentation in R/wrapper_dds.R
+\name{wrap_dds}
+\alias{wrap_dds}
 \title{Wrapper Function for DESeq2 Analysis}
 \usage{
-wrapper_DESeq2(dds, lfcThreshold, altHypothesis, correction_method = "BH")
+wrap_dds(dds, lfcThreshold, altHypothesis, correction_method = "BH")
 }
 \arguments{
 \item{dds}{A DESeqDataSet object containing the count data and experimental design.}
@@ -38,5 +38,5 @@ mock_data <- mock_rnaseq(input_var_list, N_GENES, MIN_REPLICATES, MAX_REPLICATES
 dds <- DESeq2::DESeqDataSetFromMatrix(mock_data$counts , 
                    mock_data$metadata, ~ genotype + environment)
 dds <- DESeq2::DESeq(dds, quiet = TRUE)
-result <- wrapper_DESeq2(dds, lfcThreshold = 1, altHypothesis = "greater")
+result <- wrap_dds(dds, lfcThreshold = 1, altHypothesis = "greater")
 }
diff --git a/tests/testthat/test-actualinteractionfixeffects.R b/tests/testthat/test-actual_interactionfixeffects.R
similarity index 100%
rename from tests/testthat/test-actualinteractionfixeffects.R
rename to tests/testthat/test-actual_interactionfixeffects.R
diff --git a/tests/testthat/test-actualmainfixeffects.R b/tests/testthat/test-actual_mainfixeffects.R
similarity index 93%
rename from tests/testthat/test-actualmainfixeffects.R
rename to tests/testthat/test-actual_mainfixeffects.R
index f8d8d3279749fe982684ae593cd22e36cba58b75..0c47b6aae29a7214fb4163570b124c188708ba1e 100644
--- a/tests/testthat/test-actualmainfixeffects.R
+++ b/tests/testthat/test-actual_mainfixeffects.R
@@ -47,23 +47,7 @@ test_that("averageByGroup returns correct average values", {
   expect_equal(result$logQij_mean, c(1, 3,5, 7, 2, 4, 6, 8))  # Average values
 })
 
-# Tests for convert2Factor
-test_that("convert2Factor converts specified columns to factors", {
-  # Create a sample data frame
-  data <- data.frame(
-    Category1 = c("A", "B", "A", "B"),
-    Category2 = c("X", "Y", "X", "Z"),
-    Value = 1:4,
-    stringsAsFactors = FALSE
-  )
-  
-  # Convert columns to factors
-  result <- convert2Factor(data, columns = c("Category1", "Category2"))
-  
-  # Check the output
-  expect_is(result$Category1, "factor")  # Category1 column converted to factor
-  expect_is(result$Category2, "factor")  # Category2 column converted to factor
-})
+
 
 # Tests for findAttribute
 test_that("findAttribute returns the correct attribute", {
@@ -220,4 +204,3 @@ test_that("generateActualForMainFixEff returns correct values for main fixed eff
   expect_equal(df_term, expected)
 })
 
-
diff --git a/tests/testthat/test-scalinggeneexpression.R b/tests/testthat/test-basal_expression__scaling.R
similarity index 90%
rename from tests/testthat/test-scalinggeneexpression.R
rename to tests/testthat/test-basal_expression__scaling.R
index 93f1391490532594c25b9b1cc7814fe8a77e48d3..86876f855ef01a1028d16ca8b40b76950f177c14 100644
--- a/tests/testthat/test-scalinggeneexpression.R
+++ b/tests/testthat/test-basal_expression__scaling.R
@@ -1,15 +1,15 @@
 # WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
 
 
-test_that("generate_BE returns correct number of genes", {
-  be_data <- generate_BE(n_genes = 100, 1)
+test_that("generate_basal_expression returns correct number of genes", {
+  be_data <- generate_basal_expression(n_genes = 100, 1)
   expect_equal(nrow(be_data), 100)
 })
 
 
-test_that("generate_BE returns BE values within specified vector", {
+test_that("generate_basal_expression returns BE values within specified vector", {
   BE_vec <- c(1, 2, 33, 0.4)
-  be_data <- generate_BE(n_genes = 100, BE_vec)
+  be_data <- generate_basal_expression(n_genes = 100, BE_vec)
   expect_true(all(be_data$basalExpr %in% BE_vec))
 })
 
diff --git a/tests/testthat/test-countsplot.R b/tests/testthat/test-counts_plot.R
similarity index 100%
rename from tests/testthat/test-countsplot.R
rename to tests/testthat/test-counts_plot.R
diff --git a/tests/testthat/test-evaluatedispersion.R b/tests/testthat/test-evaluate_dispersion.R
similarity index 87%
rename from tests/testthat/test-evaluatedispersion.R
rename to tests/testthat/test-evaluate_dispersion.R
index fef80adca74ccccb05a34cc2ec0af36d67699699..76e69a222ac3ec0cfb5a966ba165de73e9e613e3 100644
--- a/tests/testthat/test-evaluatedispersion.R
+++ b/tests/testthat/test-evaluate_dispersion.R
@@ -16,7 +16,7 @@ test_that("dispersion_plot function works correctly", {
   expect_s3_class(disp_plot, "gg")
 })
 
-test_that("extractTMBDispersion function extracts dispersion correctly", {
+test_that("extract_tmbDispersion function extracts dispersion correctly", {
    N_GENES = 100
   MAX_REPLICATES = 5
   MIN_REPLICATES = 5
@@ -27,11 +27,11 @@ test_that("extractTMBDispersion function extracts dispersion correctly", {
   l_res <- fitModelParallel(formula = kij ~ varA,
                           data = data2fit, group_by = "geneID",
                           family = glmmTMB::nbinom2(link = "log"), n.cores = 1)
-  extracted_disp <- extractTMBDispersion(l_res)
+  extracted_disp <- extract_tmbDispersion(l_res)
   expect_identical(colnames(extracted_disp), c("inferred_dispersion", "geneID"))
 })
 
-test_that("extractDESeqDispersion function extracts dispersion correctly", {
+test_that("extract_ddsDispersion function extracts dispersion correctly", {
    N_GENES = 100
   MAX_REPLICATES = 5
   MIN_REPLICATES = 5
@@ -43,9 +43,9 @@ test_that("extractDESeqDispersion function extracts dispersion correctly", {
       colData = mock_data$metadata,
       design = ~ varA)
   dds <- DESeq2::DESeq(dds, quiet = TRUE)
-  deseq_wrapped = wrapper_DESeq2(dds, 2, "greaterAbs")
+  deseq_wrapped = wrap_dds(dds, 2, "greaterAbs")
   
-  extracted_disp <- extractDESeqDispersion(deseq_wrapped)
+  extracted_disp <- extract_ddsDispersion(deseq_wrapped)
   expect_identical(colnames(extracted_disp), c("inferred_dispersion", "geneID"))
 })
 
@@ -61,7 +61,7 @@ test_that("getDispersionComparison function works correctly", {
                           data = data2fit, group_by = "geneID",
                           family = glmmTMB::nbinom2(link = "log"), n.cores = 1)
   
-  tmb_disp_inferred <- extractTMBDispersion(l_res)
+  tmb_disp_inferred <- extract_tmbDispersion(l_res)
     
   comparison <- getDispersionComparison(tmb_disp_inferred, mock_data$groundTruth$gene_dispersion)
   expect_identical(colnames(comparison), c("actual_dispersion",  "geneID", "inferred_dispersion"))
@@ -83,12 +83,12 @@ test_that("evaluateDispersion function works correctly", {
       colData = mock_data$metadata,
       design = ~ varA)
   dds <- DESeq2::DESeq(dds, quiet = TRUE)
-  deseq_wrapped = wrapper_DESeq2(dds, 2, "greaterAbs")
+  deseq_wrapped = wrap_dds(dds, 2, "greaterAbs")
   
-  tmb_disp_inferred <- extractTMBDispersion(l_res)
+  tmb_disp_inferred <- extract_tmbDispersion(l_res)
   TMB_dispersion_df <- getDispersionComparison(tmb_disp_inferred, mock_data$groundTruth$gene_dispersion)
   TMB_dispersion_df$from <- 'HTRfit'
-  DESEQ_disp_inferred <- extractDESeqDispersion(deseq_wrapped)
+  DESEQ_disp_inferred <- extract_ddsDispersion(deseq_wrapped)
   DESEQ_dispersion_df <- getDispersionComparison(DESEQ_disp_inferred , mock_data$groundTruth$gene_dispersion)
   DESEQ_dispersion_df$from <- 'DESeq2'
     
diff --git a/tests/testthat/test-evaluationwithmixedeffect.R b/tests/testthat/test-evaluation_withmixedeffect.R
similarity index 100%
rename from tests/testthat/test-evaluationwithmixedeffect.R
rename to tests/testthat/test-evaluation_withmixedeffect.R
diff --git a/tests/testthat/test-fitmodel.R b/tests/testthat/test-fitmodel.R
index 0ea0fb5498f8e2cf86ecd55aa86e2a0aa362b347..5307cf2a4640386986675ed18b2e901f57a63c08 100644
--- a/tests/testthat/test-fitmodel.R
+++ b/tests/testthat/test-fitmodel.R
@@ -16,28 +16,28 @@ test_that("isValidInput2fit raises an error for missing variable", {
   expect_error(isValidInput2fit(iris, formula), "Variable NonExistentVariable not found")
 })
 
-test_that(".fitModel returns a fitted model object", {
+test_that("fitModel returns a fitted model object", {
   data(mtcars)
   formula <- mpg ~ cyl + disp
-  fitted_model <- suppressWarnings(.fitModel(formula, mtcars))
-  #expect_warning(.fitModel(formula, mtcars))
+  fitted_model <- suppressWarnings(fitModel(formula, mtcars))
+  #expect_warning(fitModel(formula, mtcars))
   expect_s3_class(fitted_model, "glmmTMB")
   
   # Test with invalid formula
   invalid_formula <- mpg ~ cyl + disp + invalid_var
-  expect_error(.fitModel(invalid_formula, mtcars))
+  expect_error(fitModel(invalid_formula, mtcars))
   
   
    # Additional parameters: 
    #change family + formula
   formula <- Sepal.Length ~ Sepal.Width + Petal.Length + (1 | Species)
-  fitted_models <- suppressWarnings(.fitModel(formula = formula, 
+  fitted_models <- suppressWarnings(fitModel(formula = formula, 
                                                     data = iris, 
                                                     family = glmmTMB::nbinom1(link = "log") ))
   expect_s3_class(fitted_models$call$family, "family")
   expect_equal(fitted_models$call$formula, formula)
   #change control settings
-  fitted_models <- suppressWarnings(.fitModel(formula = formula, 
+  fitted_models <- suppressWarnings(fitModel(formula = formula, 
                                                     data = iris, 
                                                     family = glmmTMB::nbinom1(link = "log"), 
                                                 control = glmmTMB::glmmTMBControl(optCtrl=list(iter.max=1e3,
@@ -122,23 +122,23 @@ test_that("Identify full-rank model matrix (with random eff)", {
 #  expect_null(fitted_model_error)
 #})
 
-test_that(".subsetData_andfit returns a glmTMB obj", {
+test_that("subsetData_andfit returns a glmTMB obj", {
   data(iris)
   group <- "setosa"
   group_by <- "Species"
   formula <- Sepal.Length ~ Sepal.Width + Petal.Length
-  fitted_model <- .subsetData_andfit(group, group_by, formula, iris)
+  fitted_model <- subsetData_andfit(group, group_by, formula, iris)
   expect_s3_class(fitted_model, "glmmTMB")
 
   # Test with invalid formula
   invalid_formula <- Sepal.Length ~ Sepal.Width + Petal.Length +  invalid_var
-  expect_error(.subsetData_andfit(group, group_by, invalid_formula, mtcars))
+  expect_error(subsetData_andfit(group, group_by, invalid_formula, mtcars))
   
   
     # Additional parameters: 
    #change family + formula
   formula <- Sepal.Length ~ Sepal.Width + Petal.Length + (1 | Species)
-  fitted_models <- suppressWarnings(.subsetData_andfit(group,
+  fitted_models <- suppressWarnings(subsetData_andfit(group,
                                                        group_by,
                                                        formula = formula, 
                                                         data = iris, 
@@ -146,7 +146,7 @@ test_that(".subsetData_andfit returns a glmTMB obj", {
   expect_s3_class(fitted_models$call$family, "family")
   expect_equal(fitted_models$call$formula, formula)
   #change control settings
-  fitted_models <- suppressWarnings(.subsetData_andfit(group,
+  fitted_models <- suppressWarnings(subsetData_andfit(group,
                                                        group_by,
                                                        formula = formula, 
                                                         data = iris, 
@@ -192,18 +192,18 @@ test_that("launchFit handles warnings and errors during the fitting process", {
   expect_equal(fitted_models$call$control,  glmmTMB::glmmTMBControl(optCtrl=list(iter.max=1e3,eval.max=1e3)))
 })
 
-test_that(".parallel_fit returns a list of fitted model objects or NULL for any errors", {
+test_that("parallel_fit returns a list of fitted model objects or NULL for any errors", {
   data(iris)
   groups <- unique(iris$Species)
   group_by <- "Species"
   formula <- Sepal.Length ~ Sepal.Width + Petal.Length
-  fitted_models <- .parallel_fit(groups, group_by, formula, iris, log_file = "log.txt", n.cores = 1)
+  fitted_models <- parallel_fit(groups, group_by, formula, iris, log_file = "log.txt", n.cores = 1)
   expect_s3_class(fitted_models$setosa, "glmmTMB")
   expect_length(fitted_models, length(groups))
 
   # Test with invalid formula
   invalid_formula <- blabla ~ cyl + disp 
-  result <- suppressWarnings(.parallel_fit(groups, group_by, invalid_formula,  
+  result <- suppressWarnings(parallel_fit(groups, group_by, invalid_formula,  
                                            iris, log_file = "log.txt",  n.cores = 1))
   expect_equal(result, expected = list(setosa = NULL, versicolor = NULL, virginica = NULL))
   
@@ -211,7 +211,7 @@ test_that(".parallel_fit returns a list of fitted model objects or NULL for any
    # Additional parameters: 
    #change family + formula
   formula <- Sepal.Length ~ Sepal.Width + Petal.Length
-  fitted_models <- suppressWarnings(.parallel_fit(formula = formula, 
+  fitted_models <- suppressWarnings(parallel_fit(formula = formula, 
                                                     data = iris, 
                                                     group_by = group_by, 
                                                     groups = "setosa",
@@ -221,7 +221,7 @@ test_that(".parallel_fit returns a list of fitted model objects or NULL for any
   expect_s3_class(fitted_models$setosa$call$family, "family")
   expect_equal(fitted_models$setosa$call$formula, formula)
   #change control settings
-  fitted_models <- suppressWarnings(.parallel_fit(formula = formula, 
+  fitted_models <- suppressWarnings(parallel_fit(formula = formula, 
                                                     data = iris, 
                                                     group_by = group_by, 
                                                     groups = "setosa",
diff --git a/tests/testthat/test-glance_tmb.R b/tests/testthat/test-glance_glmmtmb.R
similarity index 100%
rename from tests/testthat/test-glance_tmb.R
rename to tests/testthat/test-glance_glmmtmb.R
diff --git a/tests/testthat/test-identityplot.R b/tests/testthat/test-identity_plot.R
similarity index 100%
rename from tests/testthat/test-identityplot.R
rename to tests/testthat/test-identity_plot.R
diff --git a/tests/testthat/test-mock-rnaseq.R b/tests/testthat/test-mock_rnaseq.R
similarity index 54%
rename from tests/testthat/test-mock-rnaseq.R
rename to tests/testthat/test-mock_rnaseq.R
index 5f8bba26b3390fc6d07aad3817980c3f17733de6..f08ef4cdaf38ca0ddaa442f6220eb32f2d451b24 100644
--- a/tests/testthat/test-mock-rnaseq.R
+++ b/tests/testthat/test-mock_rnaseq.R
@@ -1,57 +1,17 @@
 # WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
 
 
-# Test case: Valid input vector with numeric and positive values
-test_that("Valid input vector with numeric and positive values", {
-  input_vector <- c(0.5, 1.2, 0.8)
-  result <- getValidDispersion(input_vector)
-  expect_identical(result, input_vector)
-})
-
-# Test case: Valid input vector with numeric, positive, and negative values
-test_that("Valid input vector with numeric, positive, and negative values", {
-  input_vector <- c(0.5, -0.3, 1.2, 0.8)
-  result <- getValidDispersion(input_vector)
-  expect_identical(result, c(0.5, 1.2, 0.8))
-})
-
-# Test case: Valid input vector with non-numeric elements
-test_that("Valid input vector with non-numeric elements", {
-  input_vector <- c(0.5, "invalid", 0.8)
-  result <- getValidDispersion(input_vector)
-  expect_identical(result, c(0.5, 0.8))
-})
-
-# Test case: Empty input vector
-test_that("Empty input vector", {
-  input_vector <- numeric(0)
-  expect_error(getValidDispersion(input_vector), "Invalid dispersion values provided.")
-})
-
-# Test case: unique value in vector
-test_that("unique value in vector", {
-  input_vector <- 5
-  expect_equal(getValidDispersion(input_vector), 5)
-})
-
-# Test case: All negative values
-test_that("All negative values", {
-  input_vector <- c(-0.5, -1.2, -0.8)
-  expect_error(getValidDispersion(input_vector), "Invalid dispersion values provided.")
-})
-
-
-# Test for .isDispersionMatrixValid
-test_that(".isDispersionMatrixValid returns TRUE for valid dimensions", {
+# Test for is_dispersionMatrixValid
+test_that("is_dispersionMatrixValid returns TRUE for valid dimensions", {
   matx_dispersion <- matrix(1:6, nrow = 2, ncol = 3)
   matx_bool_replication <- matrix(TRUE, nrow = 2, ncol = 3)
-  expect_true(.isDispersionMatrixValid(matx_dispersion, matx_bool_replication))
+  expect_true(is_dispersionMatrixValid(matx_dispersion, matx_bool_replication))
 })
 
-test_that(".isDispersionMatrixValid throws an error for invalid dimensions", {
+test_that("is_dispersionMatrixValid throws an error for invalid dimensions", {
   matx_dispersion <- matrix(1:4, nrow = 2, ncol = 2)
   matx_bool_replication <- matrix(TRUE, nrow = 2, ncol = 3)
-  expect_false(.isDispersionMatrixValid(matx_dispersion, matx_bool_replication))
+  expect_false(is_dispersionMatrixValid(matx_dispersion, matx_bool_replication))
 })
 
 # Test for generateCountTable
@@ -64,16 +24,15 @@ test_that("generateCountTable generates count table with correct dimensions", {
 
 
 
-# Test for .replicateMatrix
-test_that(".replicateMatrix replicates matrix correctly", {
+# Test for replicateMatrix
+test_that("replicateMatrix replicates matrix correctly", {
   matrix <- matrix(1:9, nrow = 3, ncol = 3)
   replication_matrix <- matrix(TRUE, nrow = 3, ncol = 3)
-  replicated_matrix <- .replicateMatrix(matrix, replication_matrix)
+  replicated_matrix <- replicateMatrix(matrix, replication_matrix)
   expect_equal(dim(replicated_matrix), c(3, 9))
 })
 
 
-
 # Test for mock_rnaseq
 #test_that("mock_rnaseq returns expected output", {
   # Set up input variables
diff --git a/tests/testthat/test-prepare_data2fit.R b/tests/testthat/test-prepare_data2fit.R
index 08a7172d4905841c7fd21523e7d115fb5fa5de25..b288e616b40823007a427bef587ed4245e6c0719 100644
--- a/tests/testthat/test-prepare_data2fit.R
+++ b/tests/testthat/test-prepare_data2fit.R
@@ -28,7 +28,7 @@ test_that("getColumnWithSampleID returns column name with sampleID", {
   expected_output <- "sampleID"
   
   # Get column name with sampleID
-  column_name <- .getColumnWithSampleID(dtf_countLong, mock_data$metadata)
+  column_name <- getColumnWithSampleID(dtf_countLong, mock_data$metadata)
   
   # Check if the output matches the expected output
   expect_identical(column_name, expected_output)
diff --git a/tests/testthat/test-rocplot.R b/tests/testthat/test-roc_plot.R
similarity index 100%
rename from tests/testthat/test-rocplot.R
rename to tests/testthat/test-roc_plot.R
diff --git a/tests/testthat/test-scalingsequencingdepth.R b/tests/testthat/test-sequencing_depth_scaling.R
similarity index 100%
rename from tests/testthat/test-scalingsequencingdepth.R
rename to tests/testthat/test-sequencing_depth_scaling.R
diff --git a/tests/testthat/test-simulation.R b/tests/testthat/test-simulation.R
index 3d29ae51a9f7a4f8ac3f1b1a1634d5f9d91dd14c..c39d27f36ca4b2abea38be2317ad646c2ba5e41a 100644
--- a/tests/testthat/test-simulation.R
+++ b/tests/testthat/test-simulation.R
@@ -1,7 +1,5 @@
 # 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", {
   list_var <- init_variable()
@@ -58,3 +56,172 @@ test_that("getSubCountsTable returns the correct output", {
 })
 
 
+test_that("getReplicationMatrix returns the correct replication matrix", {
+  minN <- 2
+  maxN <- 4
+  n_samples <- 3
+  expected <- matrix(c(TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, TRUE, TRUE, TRUE, FALSE), nrow = maxN)
+  
+  set.seed(123)
+  result <- getReplicationMatrix(minN, maxN, n_samples)
+  
+  expect_equal(result, expected)
+})
+
+test_that("getSampleID return the correct list of sampleID",{
+   expect_equal(getSampleID(init_variable()), c("myVariable1", "myVariable2"))
+})
+
+# Create a test case for getMu_ij
+test_that("getMu_ij returns the correct output", {
+  # Create a sample coefficient data frame
+  dtf_coef <- data.frame(
+    log_qij = c(1, 9, 0.1),
+    basalExpr = c(2, 3, 4)
+  )
+
+    # Call the getMu_ij function
+  result <- getMu_ij(dtf_coef)
+
+  # Check if the mu_ij column is added
+  expect_true("mu_ij" %in% colnames(result))
+
+  # Check the values of mu_ij
+  #expected_mu_ij <- c(20.08554, 162754.79142 , 60.34029)
+  #expect_equal(result$mu_ij, expected_mu_ij, tolerance = 0.000001)
+})
+
+
+# Create a test case for getLog_qij
+test_that("getLog_qij returns the correct output", {
+  # Create a sample coefficient data frame
+  dtf_coef <- data.frame(
+    beta1 = c(1.2, 2.3, 3.4),
+    beta2 = c(0.5, 1.0, 1.5),
+    non_numeric = c("a", "b", "c")
+  )
+
+  # Call the getLog_qij function
+  result <- getLog_qij(dtf_coef)
+
+  # Check if the log_qij column is added
+  expect_true("log_qij" %in% colnames(result))
+
+  # Check the values of log_qij
+  expected_log_qij <- c(1.7, 3.3, 4.9)
+  expect_equal(result$log_qij, expected_log_qij)
+})
+
+test_that("getCountsTable returns the correct counts table", {
+  mat_mu_ij <- matrix(c(1,2,3,4,5,6), ncol = 3, byrow = T)
+  rownames(mat_mu_ij) <- c("gene1", "gene2")
+  colnames(mat_mu_ij) <- c("sample1", "sample2", "sample3")
+  mat_disp <- matrix(c(0.3,0.3,0.3, 0.5,0.5,0.5), ncol = 3, byrow = T)
+  rownames(mat_disp) <- c("gene1", "gene2")
+  colnames(mat_disp) <- c("sample1", "sample2", "sample3")
+  mat_repl <- matrix(c(TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE), ncol = 3, byrow = T)
+  
+  expected_df <- matrix(c(0,0,1,0,0,0,0,1,0,2,34,18,0,0,3,10,7,2), nrow = 2, byrow = T) %>% as.data.frame()
+  rownames(expected_df) <- c("gene1", "gene2")
+  colnames(expected_df) <- c("sample1_1", "sample2_1", "sample3_1", "sample1_2", 
+                             "sample2_2","sample3_2","sample1_3", "sample2_3" ,"sample3_3")
+  
+  set.seed(123)
+  result <- getCountsTable(mat_mu_ij, mat_disp, mat_repl)
+
+  expect_true(is.data.frame(result))
+  expect_equal(colnames(result), colnames(expected_df))
+  expect_equal(rownames(result), rownames(expected_df))
+
+})
+
+
+
+test_that("getSampleMetadata returns expected output", {
+  # Set up input variables
+  list_var <- init_variable()
+  n_genes <- 3
+  replicationMatrix <- matrix(TRUE, nrow = 2, ncol = 2)
+
+  # Run the function
+  result <- getSampleMetadata(list_var, n_genes, replicationMatrix)
+  
+  # Define expected output
+  expected_colnames <- c("myVariable", "sampleID")
+  expect_equal(colnames(result), expected_colnames)
+  
+  # Check the output class
+  expect_true(is.data.frame(result))
+  
+  # check nrow output
+  expect_equal(nrow(result), 4)
+
+})
+
+
+test_that("replicateByGroup return the correct ouptut", {
+  df <- data.frame(group = c("A", "B"), value = c(1, 2))
+  result <- replicateByGroup(df, "group", c(2, 3))
+  
+  expect <- data.frame(group = c("A", "A", "B", "B", "B"), 
+                       value = c(1, 1, 2,2,2), 
+                       sampleID = c("_1", "_2", "_1", "_2", "_3" ))
+  expect_equal(result, expect)
+
+})
+
+
+test_that("getDispersionMatrix returns the correct dispersion matrix", {
+  n_genes = 3
+  list_var = init_variable()
+  dispersion <- 1:3
+  expected <- matrix(1:3,byrow = F, nrow = 3, ncol = 2)
+  rownames(expected) <- c("gene1", "gene2", "gene3")
+  colnames(expected) <- c("myVariable1", "myVariable2")
+  result <- getDispersionMatrix(list_var, n_genes, dispersion )
+  expect_equal(result, expected)
+})
+
+
+# Test case: Valid input vector with numeric and positive values
+test_that("Valid input vector with numeric and positive values", {
+  input_vector <- c(0.5, 1.2, 0.8)
+  result <- getValidDispersion(input_vector)
+  expect_identical(result, input_vector)
+})
+
+# Test case: Valid input vector with numeric, positive, and negative values
+test_that("Valid input vector with numeric, positive, and negative values", {
+  input_vector <- c(0.5, -0.3, 1.2, 0.8)
+  result <- getValidDispersion(input_vector)
+  expect_identical(result, c(0.5, 1.2, 0.8))
+})
+
+# Test case: Valid input vector with non-numeric elements
+test_that("Valid input vector with non-numeric elements", {
+  input_vector <- c(0.5, "invalid", 0.8)
+  result <- getValidDispersion(input_vector)
+  expect_identical(result, c(0.5, 0.8))
+})
+
+# Test case: Empty input vector
+test_that("Empty input vector", {
+  input_vector <- numeric(0)
+  expect_error(getValidDispersion(input_vector), "Invalid dispersion values provided.")
+})
+
+# Test case: unique value in vector
+test_that("unique value in vector", {
+  input_vector <- 5
+  expect_equal(getValidDispersion(input_vector), 5)
+})
+
+# Test case: All negative values
+test_that("All negative values", {
+  input_vector <- c(-0.5, -1.2, -0.8)
+  expect_error(getValidDispersion(input_vector), "Invalid dispersion values provided.")
+})
+
+
+
+
diff --git a/tests/testthat/test-simulation2.R b/tests/testthat/test-simulation2.R
deleted file mode 100644
index 79662fbfde371c3cfc2d060f269fc9cec74c4b1e..0000000000000000000000000000000000000000
--- a/tests/testthat/test-simulation2.R
+++ /dev/null
@@ -1,131 +0,0 @@
-# WARNING - Generated by {fusen} from /dev/flat_full.Rmd: do not edit by hand
-
-
-test_that("getReplicationMatrix returns the correct replication matrix", {
-  minN <- 2
-  maxN <- 4
-  n_samples <- 3
-  expected <- matrix(c(TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, TRUE, TRUE, TRUE, FALSE), nrow = maxN)
-  
-  set.seed(123)
-  result <- getReplicationMatrix(minN, maxN, n_samples)
-  
-  expect_equal(result, expected)
-})
-
-test_that("getSampleID return the correct list of sampleID",{
-   expect_equal(getSampleID(init_variable()), c("myVariable1", "myVariable2"))
-})
-
-# Create a test case for getMu_ij
-test_that("getMu_ij returns the correct output", {
-  # Create a sample coefficient data frame
-  dtf_coef <- data.frame(
-    log_qij = c(1, 9, 0.1),
-    basalExpr = c(2, 3, 4)
-  )
-
-    # Call the getMu_ij function
-  result <- getMu_ij(dtf_coef)
-
-  # Check if the mu_ij column is added
-  expect_true("mu_ij" %in% colnames(result))
-
-  # Check the values of mu_ij
-  #expected_mu_ij <- c(20.08554, 162754.79142 , 60.34029)
-  #expect_equal(result$mu_ij, expected_mu_ij, tolerance = 0.000001)
-})
-
-
-# Create a test case for getLog_qij
-test_that("getLog_qij returns the correct output", {
-  # Create a sample coefficient data frame
-  dtf_coef <- data.frame(
-    beta1 = c(1.2, 2.3, 3.4),
-    beta2 = c(0.5, 1.0, 1.5),
-    non_numeric = c("a", "b", "c")
-  )
-
-  # Call the getLog_qij function
-  result <- getLog_qij(dtf_coef)
-
-  # Check if the log_qij column is added
-  expect_true("log_qij" %in% colnames(result))
-
-  # Check the values of log_qij
-  expected_log_qij <- c(1.7, 3.3, 4.9)
-  expect_equal(result$log_qij, expected_log_qij)
-})
-
-test_that("getCountsTable returns the correct counts table", {
-  mat_mu_ij <- matrix(c(1,2,3,4,5,6), ncol = 3, byrow = T)
-  rownames(mat_mu_ij) <- c("gene1", "gene2")
-  colnames(mat_mu_ij) <- c("sample1", "sample2", "sample3")
-  mat_disp <- matrix(c(0.3,0.3,0.3, 0.5,0.5,0.5), ncol = 3, byrow = T)
-  rownames(mat_disp) <- c("gene1", "gene2")
-  colnames(mat_disp) <- c("sample1", "sample2", "sample3")
-  mat_repl <- matrix(c(TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE), ncol = 3, byrow = T)
-  
-  expected_df <- matrix(c(0,0,1,0,0,0,0,1,0,2,34,18,0,0,3,10,7,2), nrow = 2, byrow = T) %>% as.data.frame()
-  rownames(expected_df) <- c("gene1", "gene2")
-  colnames(expected_df) <- c("sample1_1", "sample2_1", "sample3_1", "sample1_2", 
-                             "sample2_2","sample3_2","sample1_3", "sample2_3" ,"sample3_3")
-  
-  set.seed(123)
-  result <- getCountsTable(mat_mu_ij, mat_disp, mat_repl)
-
-  expect_true(is.data.frame(result))
-  expect_equal(colnames(result), colnames(expected_df))
-  expect_equal(rownames(result), rownames(expected_df))
-
-})
-
-
-
-test_that("getSampleMetadata returns expected output", {
-  # Set up input variables
-  list_var <- init_variable()
-  n_genes <- 3
-  replicationMatrix <- matrix(TRUE, nrow = 2, ncol = 2)
-
-  # Run the function
-  result <- getSampleMetadata(list_var, n_genes, replicationMatrix)
-  
-  # Define expected output
-  expected_colnames <- c("myVariable", "sampleID")
-  expect_equal(colnames(result), expected_colnames)
-  
-  # Check the output class
-  expect_true(is.data.frame(result))
-  
-  # check nrow output
-  expect_equal(nrow(result), 4)
-
-})
-
-
-test_that(".replicateByGroup return the correct ouptut", {
-  df <- data.frame(group = c("A", "B"), value = c(1, 2))
-  result <- .replicateByGroup(df, "group", c(2, 3))
-  
-  expect <- data.frame(group = c("A", "A", "B", "B", "B"), 
-                       value = c(1, 1, 2,2,2), 
-                       sampleID = c("_1", "_2", "_1", "_2", "_3" ))
-  expect_equal(result, expect)
-
-})
-
-
-test_that("getDispersionMatrix returns the correct dispersion matrix", {
-  n_genes = 3
-  list_var = init_variable()
-  dispersion <- 1:3
-  expected <- matrix(1:3,byrow = F, nrow = 3, ncol = 2)
-  rownames(expected) <- c("gene1", "gene2", "gene3")
-  colnames(expected) <- c("myVariable1", "myVariable2")
-  result <- getDispersionMatrix(list_var, n_genes, dispersion )
-  expect_equal(result, expected)
-})
-
-
-
diff --git a/tests/testthat/test-simulation_initialization.R b/tests/testthat/test-simulation_initialization.R
index b0c4b1116f1fe43ad1ec68944abab47c137e34de..01dffa29db1e9084801af2016669e9bb24dedf14 100644
--- a/tests/testthat/test-simulation_initialization.R
+++ b/tests/testthat/test-simulation_initialization.R
@@ -90,13 +90,6 @@ test_that("build_sub_obj_return_to_user returns the expected output", {
   
 })
 
-test_that("generateGridCombination_fromListVar returns expected output", {
-  result <- generateGridCombination_fromListVar(init_variable())
-  expect <- data.frame(label_myVariable = c("myVariable1", "myVariable2"))
-  expect_equal(nrow(result), nrow(expect))
-  expect_equal(ncol(result), ncol(expect))
-})
-
 test_that("add_interaction adds an interaction between variables", {
   list_var <- init_variable(name = "varA", mu = 1, sd = 1, level = 2)
   list_var <- init_variable(list_var, name = "varB", mu = 2, sd = 1, level = 3)
@@ -115,17 +108,3 @@ test_that("getNumberOfCombinationsInInteraction calculates the number of combina
   list_var <- init_variable(list_var, name = "varB", mu = 2, sd = 1, level = 3)
   expect_equal(getNumberOfCombinationsInInteraction(list_var, c("varA", "varB")), 6)
 })
-
-test_that("getLabels generates labels for variables", {
-  labels <- getLabels(c("varA", "varB"), c(2, 3))
-  expect_equal(length(labels), 2)
-  expect_equal(length(labels[[1]]), 2)
-  expect_equal(length(labels[[2]]), 3)
-})
-
-test_that("getGridCombination generates a grid of combinations", {
-  labels <- list(A = c("A1", "A2"), B = c("B1", "B2", "B3"))
-  grid_combination <- getGridCombination(labels)
-  expect_equal(dim(grid_combination), c(6, 2))
-})
-
diff --git a/tests/testthat/test-simulationreport.R b/tests/testthat/test-simulation_report.R
similarity index 100%
rename from tests/testthat/test-simulationreport.R
rename to tests/testthat/test-simulation_report.R
diff --git a/tests/testthat/test-tidy_glmmtmb.R b/tests/testthat/test-tidy_glmmtmb.R
index 2b7a36d9bc2e129aba20eacd9e043b5d032b11d2..06ce9c42ca0e6b57d713a17c1e57da8f9df6283b 100644
--- a/tests/testthat/test-tidy_glmmtmb.R
+++ b/tests/testthat/test-tidy_glmmtmb.R
@@ -53,13 +53,6 @@ test_that("build_missingColumn_with_na returns the correct results", {
 })
 
 
-test_that("removeDuplicatedWord returns expected output", {
-  words <- c("hellohello", "worldworld", "programmingprogramming", "R isis great")
-  cleaned_words <- removeDuplicatedWord(words)
-  expect_equal(cleaned_words, c("hello", "world", "programming", "R is great"))
-})
-
-
 
 test_that("correlation_matrix_2df returns expected output",{
 
@@ -107,19 +100,6 @@ test_that("renameColumns returns expected output",{
   expect_equal(dim(renamed_df), dim(df))
 })
     
-
-test_that("reorderColumns returns expected output",{
-    df <- data.frame(A = 1:3, B = 4:6, C = 7:9)
-    # Define the desired column order
-    columnOrder <- c("B", "C", "A")
-    # Reorder the columns of the dataframe
-    df_reorder <- reorderColumns(df, columnOrder)
-    expect_equal(colnames(df_reorder), columnOrder)
-    expect_equal(dim(df_reorder), dim(df))
-
-})
-
-
 test_that("tidy_tmb returns expected output",{
   model1 <- glmmTMB::glmmTMB(Sepal.Length ~ Sepal.Width + Petal.Length + (1 | Species), data = iris)
   model2 <- glmmTMB::glmmTMB(Petal.Length ~ Sepal.Length + Sepal.Width + (1 | Species), data = iris)
diff --git a/tests/testthat/test-updatefitmodel.R b/tests/testthat/test-update_fittedmodel.R
similarity index 90%
rename from tests/testthat/test-updatefitmodel.R
rename to tests/testthat/test-update_fittedmodel.R
index b7f15d5305acac42d061390f76b3bdc3e73c7e0d..f38b3080449cf266d037ee57c4e3e42b330bd454 100644
--- a/tests/testthat/test-updatefitmodel.R
+++ b/tests/testthat/test-update_fittedmodel.R
@@ -21,14 +21,14 @@ test_that("updateParallel function returns correct results", {
   # Additional parameters: 
    #change family + formula
   new_formula <- Sepal.Length ~ Sepal.Width 
-  updated_model <- suppressWarnings(updateParallel(l_tmb = fitted_models, 
+  updated_model <- suppressWarnings(updateParallel(fitted_models, 
                                                     formula = new_formula,
                                                     n.cores = 1,
                                                     family = glmmTMB::nbinom1(link = "log") ))
   expect_s3_class(updated_model$setosa$call$family, "family")
   expect_equal(updated_model$setosa$call$formula, new_formula)
   #change control settings
-  updated_model <- suppressWarnings(updateParallel(l_tmb = fitted_models, 
+  updated_model <- suppressWarnings(updateParallel(fitted_models, 
                                                  formula = new_formula, 
                                                  family = glmmTMB::nbinom1(link = "log"), 
                                                   n.cores = 1,
@@ -37,15 +37,15 @@ test_that("updateParallel function returns correct results", {
   expect_equal(updated_model$setosa$call$control,  glmmTMB::glmmTMBControl(optCtrl=list(iter.max=1e3,eval.max=1e3)))
   
   # Update an updated model
-  updated_updated_model <- suppressWarnings(updateParallel(l_tmb = updated_model, 
+  updated_updated_model <- suppressWarnings(updateParallel(updated_model, 
                                                  formula = new_formula, 
                                                   n.cores = 1,
                                                  family = glmmTMB::ziGamma(link = "inverse")))
   expect_s3_class(updated_updated_model$setosa$call$family,  "family")
 })
 
-# Test .parallel_update function
-test_that(".parallel_update function returns correct results", {
+# Test parallel_update function
+test_that("parallel_update function returns correct results", {
 # Load the required data
   data(iris)
   groups <- unique(iris$Species)
@@ -53,14 +53,14 @@ test_that(".parallel_update function returns correct results", {
   formula <- Sepal.Length ~ Sepal.Width + Petal.Length
   fitted_models <- fitModelParallel(formula, iris, group_by, n.cores = 1)
   new_formula <- Sepal.Length ~ Sepal.Width 
-  results <- .parallel_update(new_formula, fitted_models, n.cores = 1)
+  results <- parallel_update(new_formula, fitted_models, n.cores = 1)
   expect_is(results, "list")
   expect_equal(length(results), length(fitted_models))
   expect_is(results$setosa, "glmmTMB")
 
   #uncorrect formula 
   new_formula <- Sepal.Length ~ blabla
-  results <- .parallel_update(new_formula, fitted_models, n.cores = 1)
+  results <- parallel_update(new_formula, fitted_models, n.cores = 1)
   expect_is(results, "list")
   expect_equal(length(results), length(fitted_models))
   expect_equal(results$setosa, NULL)
@@ -68,14 +68,14 @@ test_that(".parallel_update function returns correct results", {
   # Additional parameters: 
    #change family + formula
   new_formula <- Sepal.Length ~ Sepal.Width 
-  updated_model <- suppressWarnings(.parallel_update(l_tmb = fitted_models, 
+  updated_model <- suppressWarnings(parallel_update(fitted_models, 
                                                      formula = new_formula,
                                                       n.cores = 1,
                                                       family = glmmTMB::nbinom1(link = "log") ))
   expect_s3_class(updated_model$setosa$call$family, "family")
   expect_equal(updated_model$setosa$call$formula, new_formula)
   #change control
-  updated_model <- suppressWarnings(.parallel_update(l_tmb = fitted_models, 
+  updated_model <- suppressWarnings(parallel_update(fitted_models, 
                                                  formula = new_formula, 
                                                   n.cores = 1,
                                                  family = glmmTMB::nbinom1(link = "log"), 
diff --git a/tests/testthat/test-utils.R b/tests/testthat/test-utils.R
index 9aaf15c6080753e4c13fd037be810139847617f6..bca87e27de380be42b17f0504db736677f0d0d4b 100644
--- a/tests/testthat/test-utils.R
+++ b/tests/testthat/test-utils.R
@@ -29,3 +29,87 @@ test_that("clean_variable_name handles reserved names properly", {
   expect_error(clean_variable_name("interactions"))
   expect_error(clean_variable_name("correlations"))
 })
+
+
+test_that("getLabels generates labels for variables", {
+  labels <- getLabels(c("varA", "varB"), c(2, 3))
+  expect_equal(length(labels), 2)
+  expect_equal(length(labels[[1]]), 2)
+  expect_equal(length(labels[[2]]), 3)
+})
+
+test_that("getGridCombination generates a grid of combinations", {
+  labels <- list(A = c("A1", "A2"), B = c("B1", "B2", "B3"))
+  grid_combination <- getGridCombination(labels)
+  expect_equal(dim(grid_combination), c(6, 2))
+})
+
+
+test_that("generateGridCombination_fromListVar returns expected output", {
+  result <- generateGridCombination_fromListVar(init_variable())
+  expect <- data.frame(label_myVariable = c("myVariable1", "myVariable2"))
+  expect_equal(nrow(result), nrow(expect))
+  expect_equal(ncol(result), ncol(expect))
+})
+
+# Tests for convert2Factor
+test_that("convert2Factor converts specified columns to factors", {
+  # Create a sample data frame
+  data <- data.frame(
+    Category1 = c("A", "B", "A", "B"),
+    Category2 = c("X", "Y", "X", "Z"),
+    Value = 1:4,
+    stringsAsFactors = FALSE
+  )
+  
+  # Convert columns to factors
+  result <- convert2Factor(data, columns = c("Category1", "Category2"))
+  
+  # Check the output
+  expect_is(result$Category1, "factor")  # Category1 column converted to factor
+  expect_is(result$Category2, "factor")  # Category2 column converted to factor
+})
+
+test_that("removeDuplicatedWord returns expected output", {
+  words <- c("hellohello", "worldworld", "programmingprogramming", "R isis great")
+  cleaned_words <- removeDuplicatedWord(words)
+  expect_equal(cleaned_words, c("hello", "world", "programming", "R is great"))
+})
+
+test_that("reorderColumns returns expected output",{
+    df <- data.frame(A = 1:3, B = 4:6, C = 7:9)
+    # Define the desired column order
+    columnOrder <- c("B", "C", "A")
+    # Reorder the columns of the dataframe
+    df_reorder <- reorderColumns(df, columnOrder)
+    expect_equal(colnames(df_reorder), columnOrder)
+    expect_equal(dim(df_reorder), dim(df))
+
+})
+
+
+test_that( "generateGridCombination_fromListVar return expected output", {
+    ## case 1
+    gridcom <- generateGridCombination_fromListVar(init_variable())
+    expect_s3_class(gridcom, "data.frame")
+    expect_equal(gridcom$label_myVariable, factor(c("myVariable1", "myVariable2")))
+
+    ## case 2
+    init_variables <- init_variable() %>% init_variable(name = "var" , mu = 2, sd = 1, level = 3) 
+    gridcom <- generateGridCombination_fromListVar(init_variables)
+    expect_s3_class(gridcom, "data.frame")
+    expect_equal(unique(gridcom$label_myVariable), factor(c("myVariable1", "myVariable2")))
+    expect_equal(unique(gridcom$label_var), factor(c("var1", "var2", "var3")))
+
+  })
+
+test_that( "getGivenAttribute return expected output", {
+  ## -- case 1
+  level_attr <- getGivenAttribute(init_variable(), "level")
+  expect_equal(level_attr$myVariable, 2)
+
+  ## -- case 2
+  init_variables <- init_variable() %>% init_variable(name = "var" , mu = 2, sd = 1, level = 3) 
+  mu_attr <- getGivenAttribute(init_variables, "mu")
+  expect_equal(mu_attr$var, 2)
+} )
diff --git a/tests/testthat/test-wrapperdeseq2.R b/tests/testthat/test-wrapper_dds.R
similarity index 94%
rename from tests/testthat/test-wrapperdeseq2.R
rename to tests/testthat/test-wrapper_dds.R
index cee2492d4be8daa60c577b944701c9afac69c911..03d57e1e5009191dd0ba0dec6a6f5734020204af 100644
--- a/tests/testthat/test-wrapperdeseq2.R
+++ b/tests/testthat/test-wrapper_dds.R
@@ -2,7 +2,7 @@
 
 
 
-test_that("get_inference returns a data frame with correct columns", {
+test_that("get_inference_dds returns a data frame with correct columns", {
   # Create a sample dds_full data frame
   N_GENES = 100
   MAX_REPLICATES = 5
@@ -17,7 +17,7 @@ test_that("get_inference returns a data frame with correct columns", {
   dds_full <- S4Vectors::mcols(dds) %>% as.data.frame()
   
   # Call the function
-  inference_results <- get_inference(dds_full, lfcThreshold = 0.5, altHypothesis = "greater", correction_method = "BH")
+  inference_results <- get_inference_dds(dds_full, lfcThreshold = 0.5, altHypothesis = "greater", correction_method = "BH")
   
   # Check if the returned object is a data frame
   expect_true(is.data.frame(inference_results))
@@ -104,7 +104,7 @@ test_that("wrapperDESeq2 function works correctly", {
   mock_data <- mock_rnaseq(input_var_list, N_GENES, MIN_REPLICATES, max_replicates = MAX_REPLICATES)
   dds <- DESeq2::DESeqDataSetFromMatrix(mock_data$counts , mock_data$metadata, ~ genotype + environment)
   dds <- DESeq2::DESeq(dds, quiet = TRUE)
-  deseq2_wrapped <- wrapper_DESeq2(dds, 0.2, "greaterAbs")
+  deseq2_wrapped <- wrap_dds(dds, 0.2, "greaterAbs")
   
   expect_true(is.list(deseq2_wrapped))
   
diff --git a/vignettes/htrfit.Rmd b/vignettes/htrfit.Rmd
index 56f5478f12a7c8efb5ee54ca613fd9e990dfcec4..b2a91288d232c3a199f9a4eefe0653fc8141dc5c 100644
--- a/vignettes/htrfit.Rmd
+++ b/vignettes/htrfit.Rmd
@@ -15,7 +15,6 @@ knitr::opts_chunk$set(
 ```
 
 ```{r setup}
-devtools::load_all()
 library(HTRfit)
 ```
 
@@ -316,14 +315,14 @@ The `updateParallel()` function updates and re-fits a model for each gene. It of
 ```{r example-update, warning = FALSE, message = FALSE}
 ## -- update your fit modifying the model family
 l_tmb <- updateParallel(formula =  kij ~ varA,
-                          l_tmb = l_tmb ,
+                          list_tmb = l_tmb ,
                           family = gaussian(), 
                           log_file = "log.txt",
                           n.cores = 1)
 
 ## -- update fit using additional model control settings
 l_tmb <- updateParallel(formula =  kij ~ varA ,
-                          l_tmb = l_tmb ,
+                          list_tmb = l_tmb ,
                           family = gaussian(), 
                           log_file = "log.txt",
                           n.cores = 1,
@@ -333,7 +332,7 @@ l_tmb <- updateParallel(formula =  kij ~ varA ,
 
 ## -- update your model formula and your family model
 l_tmb <- updateParallel(formula =   kij ~ varA + varB  + varA:varB ,
-                          l_tmb = l_tmb ,
+                          list_tmb = l_tmb ,
                           family = glmmTMB::nbinom2(link = "log"), 
                           log_file = "log.txt",
                           n.cores = 1)
@@ -349,12 +348,12 @@ Visualizing fit metrics is essential for evaluating your models. Here, we show y
 
 ```{r example-plotMetrics, warning = FALSE, message = FALSE, fig.align = 'center', fig.height = 4, fig.width = 6}
 ## -- plot all metrics
-metrics_plot(l_tmb = l_tmb)
+metrics_plot(list_tmb = l_tmb)
 ```
 
 ```{r example-plotMetricsFocus, warning = FALSE, message = FALSE, fig.align = 'center', fig.height = 3, fig.width = 4}
 ## -- Focus on metrics
-metrics_plot(l_tmb = l_tmb, focus = c("dispersion", "logLik"))
+metrics_plot(list_tmb = l_tmb, focus = c("dispersion", "logLik"))
 ```
 
 ## Anova to select the best model
@@ -364,10 +363,10 @@ Utilizing the `anovaParallel()` function enables you to perform model selection
 
 ```{r example-anova, warning = FALSE, message = FALSE}
 ## -- update your fit modifying the model family
-l_anova <- anovaParallel(l_tmb = l_tmb)
+l_anova <- anovaParallel(list_tmb = l_tmb)
 
 ## -- additional settings
-l_anova <- anovaParallel(l_tmb = l_tmb, type = "III" )
+l_anova <- anovaParallel(list_tmb = l_tmb, type = "III" )
 
 ## -- output 
 l_anova$gene1