Skip to content
Snippets Groups Projects
Commit 7eb2c178 authored by Arnaud Duvermy's avatar Arnaud Duvermy
Browse files

update src file

parent 2611e4e9
No related branches found
No related tags found
No related merge requests found
############################# PCKGE REQUIRED ##############################
library(DESeq2)
library(ggplot2)
library(tydiverse)
### maybe others ###
# fix seed
set.seed(123)
......@@ -6,50 +13,17 @@ set.seed(123)
source("mydatalocal/counts_simulation/src/simulators.R")
# visualization functions
mu_effect_visualization <- function(mu_effect_res){
label_wrap <- c("mu observed", "N gene DE", "min(|logFC|)", "var observed")
names(label_wrap) <- c("mu_observ", "res_DEA", "statistical_power", "var_observ")
figure = mu_effect_res %>% ggplot(., aes(x=vec_of_mu, y = value, col=factor(N_rep))) +
geom_point() + facet_wrap(~variable, scales = "free_y", labeller = labeller(variable = label_wrap)) + labs(color = "N replicates")
return(figure)
}
size_effect_visualization <- function(alpha_effect_res){
label_wrap <- c("mu observed", "N gene DE", "min(|logFC|)", "var observed")
names(label_wrap) <- c("mu_observ", "res_DEA", "statistical_power", "var_observ")
figure = alpha_effect_res %>% ggplot(., aes(x=vec_of_alpha, y = value, col=factor(N_rep))) +
geom_point() + facet_wrap(~variable, scales = "free_y", labeller = labeller(variable = label_wrap)) + labs(color = "N replicates")
return(figure)
}
#visualization function
source("mydatalocal/counts_simulation/src/visualization_fun.R")
## main
# Params to specify design
min_rep = 2 #/!\ = 1 forbidden
max_rep = 10
########################## INPUT PARAMS #####################################
N_cond = 2
N_gene = 6000
n_rep_sim = seq(2, 5, by = 1) ### number of replicate to assessed
mu_simul = seq(100, 15000, by = 200)
mu_simul
mu_simul <- rep.int(1500, 8)
res_simul <- mu_effect(alpha = 2.9, mu_simul)
reshape_res_simul <- res_simul %>% reshape2::melt(.,id = c("vec_of_mu"))
mu_effect_visualization(reshape_res_simul)
alpha_simul = seq(0.01, 10, by = 0.1)
alpha_simul
res_simul2 <- size_effect(mu = 10000, alpha_simul)
res_simul2
reshape_res_simul2 <- res_simul2 %>% reshape2::melt(.,id = c("vec_of_alpha"))
size_effect_visualization(reshape_res_simul2)
## replicate effect
n_rep_sim = seq(2, 5, by = 1)
############################ MU effect #######################################
mu_simul_dtf_res <- data.frame()
for (N_rep in n_rep_sim){
......@@ -62,16 +36,17 @@ for (N_rep in n_rep_sim){
mu_simul_dtf_res <- rbind(mu_simul_dtf_res, tmp_reshape_res_simul)
}
## LOG transform
###### LOG transform #######
# -> SEE linearity of var observed & mu
#mu_simul_dtf_res$value[mu_simul_dtf_res$variable=="var_observ"]<-log(mu_simul_dtf_res$value[mu_simul_dtf_res$variable=="var_observ"])
#mu_simul_dtf_res$vec_of_mu <- log(mu_simul_dtf_res$vec_of_mu)
## Visualization
###### Visualization ######
figure_mu_effect <- mu_effect_visualization(mu_simul_dtf_res)
figure_mu_effect
svg("mydatalocal/counts_simulation/img/fig_mu_effect.svg")
figure_mu_effect
dev.off()
########################### ALPHA effect ####################################
n_rep_sim = seq(2, 5, by = 1)
......@@ -86,10 +61,47 @@ for (N_rep in n_rep_sim){
tmp_reshape_res_simul <- res_simul %>% reshape2::melt(.,id = c("vec_of_alpha", "N_rep"))
alpha_simul_dtf_res <- rbind(alpha_simul_dtf_res, tmp_reshape_res_simul)
}
###### Visualization ######
alpha_simul_dtf_res
figure_alpha_effect <- size_effect_visualization(alpha_simul_dtf_res)
figure_alpha_effect
########################### EXPORT RESULTS #################################
svg("mydatalocal/counts_simulation/img/fig_mu_effect.svg")
figure_mu_effect
dev.off()
svg("mydatalocal/counts_simulation/img/fig_size_effect.svg")
figure_alpha_effect
dev.off()
########################### Beta test #####################################
## main
# Params to specify design
min_rep = 2 #/!\ = 1 forbidden
max_rep = 10
mu_simul = seq(100, 15000, by = 200)
mu_simul
mu_simul <- rep.int(1500, 8)
res_simul <- mu_effect(alpha = 2.9, mu_simul)
reshape_res_simul <- res_simul %>% reshape2::melt(.,id = c("vec_of_mu"))
mu_effect_visualization(reshape_res_simul)
alpha_simul = seq(0.01, 10, by = 0.1)
alpha_simul
res_simul2 <- size_effect(mu = 10000, alpha_simul)
res_simul2
reshape_res_simul2 <- res_simul2 %>% reshape2::melt(.,id = c("vec_of_alpha"))
size_effect_visualization(reshape_res_simul2)
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment