rearrange for convenience - v3

src/v3/HTRsim/R/countsGenerator.R

@@ -215,6 +215,7 @@ uniform_replication <- function(maxN, n_samples) {
 #'
 #' @param maxN  a integer : number of replicates
 #' @param n_samples an integer : number of samples
+#' @import purrr
 #' @return a matrix of 0 and 1
 #' @export
 #'

src/v3/HTRsim/R/loadEmbeddedFiles.R

+#' Load beta dtf embedded in package
+#'
+#' @return an object dds.Extraction 
+#' @export
+#'
+#' @examples
+loadObservedValues <- function(){
+ ## Import public beta observed
+  fn = system.file("extdata/", "SRP217588_YM_observedParams.rds", package = "HTRsim")
+  dds.extraction = readRDS(file = fn)
+  return(dds.extraction)
+}
+
+#' Load public counts table
+#'
+#' @import dplyr
+#' @import utils
+#' @return dataframe
+#' @export
+#'
+#' @examples
+loadPubliCounTable <- function(){
+  ## Import public counts table
+  fn = system.file("extdata/", "SRP217588_YM_vkallisto.tsv", package = "HTRsim")
+  tabl_cnts <- utils::read.table(file = fn, header = TRUE)
+  rownames(tabl_cnts) <- tabl_cnts$gene_id
+  tabl_cnts <- tabl_cnts %>% dplyr::select(-gene_id) ##suppr colonne GeneID
+  tabl_cnts = tabl_cnts[ order(tabl_cnts %>% rownames()), ]
+  tabl_cnts = tabl_cnts %>% dplyr::select(., !matches( "ru_rm_5"))
+  return(tabl_cnts)
+}
+
+#' Load public design
+#'
+#' @import dplyr
+#' @import utils
+#' @return dataframe
+#' @export
+#'
+#' @examples
+loadPublicDesign <- function(){
+  ## DESIGN
+  fn = system.file("extdata/", "SRP217588_YM_bioDesign.csv", package = "HTRsim")
+  bioDesign <- utils::read.table(file = fn, header = T, sep = ';')
+  ## defining reference
+  bioDesign$genotype <- factor(x = bioDesign$genotype,levels = c("GSY147", "RM11"))
+  bioDesign$environment <- factor(x = bioDesign$environment, levels = c( "untreated", "treated"))
+  bioDesign = bioDesign %>% dplyr::filter(!str_detect(sample, "ru_rm_5") ) 
+  return(bioDesign)
+}
\ No newline at end of file

src/v3/HTRsim/R/manipulationsDDS_obj.R

@@ -42,12 +42,7 @@ ddsExtraction <-  function(dds_obj){
 
 
   return(list(beta = beta.dtf,
-              gene_dispersion = gene_disp,
-              beta0.mean = mean(beta.dtf$beta0,  na.rm = T), beta0.sd = sd(beta.dtf$beta0,  na.rm = T) ,
-              betaG.mean = mean(beta.dtf$betaG,  na.rm = T), betaG.sd = sd(beta.dtf$betaG,  na.rm = T) ,
-              betaE.mean = mean(beta.dtf$betaE, na.rm = T)  , betaE.sd = sd(beta.dtf$betaE,  na.rm = T) ,
-              betaGE.mean =  mean(beta.dtf$betaGE, na.rm = T)  , betaGE.sd = sd(beta.dtf$betaGE,  na.rm = T),
-              gene_disp.mean = mean(gene_disp, na.rm = T) ,  gene_disp.sd = sd(gene_disp, na.rm = T)))
+              gene_dispersion = gene_disp))
 
 }
 
@@ -60,25 +55,10 @@ ddsExtraction <-  function(dds_obj){
 #'
 #' @examples
 publicDDS_extraction <- function(){
-  ## Import & reshape public counts table
-  fn = system.file("extdata/", "SRP217588_YM_vkallisto.tsv", package = "HTRsim")
-  tabl_cnts <- utils::read.table(file = fn, header = TRUE)
-  rownames(tabl_cnts) <- tabl_cnts$gene_id
-  tabl_cnts <- tabl_cnts %>% dplyr::select(-gene_id)##suppr colonne GeneID
-  tabl_cnts = tabl_cnts[ order(tabl_cnts %>% rownames()), ]
-  ## DESIGN
-  fn = system.file("extdata/", "SRP217588_YM_bioDesign.csv", package = "HTRsim")
-  bioDesign <- utils::read.table(file = fn, header = T, sep = ';')
-
-  ## defining reference
-  bioDesign$genotype <- factor(x = bioDesign$genotype,levels = c("GSY147", "RM11"))
-  bioDesign$environment <- factor(x = bioDesign$environment, levels = c( "untreated", "treated"))
-  bioDesign = bioDesign %>% dplyr::filter(!str_detect(sample, "ru_rm_5") ) 
-  tabl_cnts = tabl_cnts %>% dplyr::select(., !matches( "ru_rm_5"))
-
+  tabl_cnts = loadPubliCounTable()
+  bioDesgin = loadPublicDesign()
   ## Launch DESEQ2
   dds = run.deseq(tabl_cnts, bioDesign = bioDesign)
-
   ## Extract
   dds.extraction = ddsExtraction(dds_obj = dds)
 

src/v3/HTRsim/R/workflow.R

@@ -4,21 +4,24 @@
 #' @param n_genes number of genes
 #' @param n_genotypes number of genotypes 
 #' @param n_environments number of env
-#' @param fit.mvnorm output of Rfast::mvnorm.mle
+#' @param dds.extraction output of dds.extraction
 #' @param max_n_replicates  max number of replicates
 #' @param uniformNumberOfReplicates boolean
 #' @param uniformDispersion boolean
-#' @return counts table 
+#' @return list 
 #' @export
 #'
 #' @examples
-counTableBuilder <- function(fit.mvnorm = getPublicMvnormFit(),
-                             n_genes, 
-                             n_genotypes, 
-                             n_environments = 2, 
-                             max_n_replicates,  
-                             uniformNumberOfReplicates = T, 
-                             uniformDispersion = T ) {
+rnaMock <- function(n_genes, 
+                    n_genotypes, 
+                    n_environments = 2, 
+                    max_n_replicates,  
+                    uniformNumberOfReplicates = T, 
+                    uniformDispersion = T,
+                    dds.extraction = loadObservedValues() ) {
+    
+    ## Fit mvnorm ##
+    fit.mvnorm  = mvnormFitting(dds.extraction$beta)
 
     ##### Ground truth ######
     beta.actual = getBetaforSimulation( n_genes, 
@@ -48,7 +51,8 @@ counTableBuilder <- function(fit.mvnorm = getPublicMvnormFit(),
                                 sample_id_list = sample_ids,  
                                 replication.matx = designReplication.matx)
     
-    return( countTable )
+    return( list( countTable = countTable, dispersion = dispersion.matrix, 
+                    beta = beta.actual, mvnorm = fit.mvnorm ))
 
 
 }