diff --git a/src/v3/HTRsim/R/countsGenerator.R b/src/v3/HTRsim/R/countsGenerator.R index a44e31e9e098a4c26a3adf9a5397a02bae57d734..6f723f6fe73756ffe798511228f5d2cbc7c44641 100644 --- a/src/v3/HTRsim/R/countsGenerator.R +++ b/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 #' diff --git a/src/v3/HTRsim/R/loadEmbeddedFiles.R b/src/v3/HTRsim/R/loadEmbeddedFiles.R new file mode 100644 index 0000000000000000000000000000000000000000..2ba458dcfc33bfda49c3fa9c2bbe57e13f01d224 --- /dev/null +++ b/src/v3/HTRsim/R/loadEmbeddedFiles.R @@ -0,0 +1,50 @@ +#' 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 diff --git a/src/v3/HTRsim/R/manipulationsDDS_obj.R b/src/v3/HTRsim/R/manipulationsDDS_obj.R index d36fa152e0839ce71ce0bccb81557460949b53ce..3689f0faf83611321a477df722fe406c32e65d79 100644 --- a/src/v3/HTRsim/R/manipulationsDDS_obj.R +++ b/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) diff --git a/src/v3/HTRsim/R/workflow.R b/src/v3/HTRsim/R/workflow.R index 41f52db0dca30024f6cba4bc099a4ab60e60f8f7..d09d2e0ab3002f03e15ad69022764916a50a9695 100644 --- a/src/v3/HTRsim/R/workflow.R +++ b/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 )) }