From 414a4f214ed1bf87dfb7fabdf6131f94d0911007 Mon Sep 17 00:00:00 2001
From: Francesco Sabatini <francesco.sabatini@idiv.de>
Date: Tue, 28 Jul 2020 16:16:28 +0200
Subject: [PATCH] Commited 02_Figures_Simulations.R

---
 02_Figures_Simulations.R | 258 +++++++++++++++++++++++++++++++++++++++
 1 file changed, 258 insertions(+)
 create mode 100644 02_Figures_Simulations.R

diff --git a/02_Figures_Simulations.R b/02_Figures_Simulations.R
new file mode 100644
index 0000000..15e9c74
--- /dev/null
+++ b/02_Figures_Simulations.R
@@ -0,0 +1,258 @@
+## this script extracts all the summary.txt files from all subfolder
+## and summarizes the output for each run, trait x environment combination, and statistics
+## It then plots the summarized output
+library(tidyverse)
+
+#### import function
+FormatData <- function(myfiles){
+  output <- NULL
+  for(ff in myfiles){
+    iter <- gsub(pattern="/Summary.txt$", replacement="", ff) 
+    iter <- strsplit(iter, split = "_")[[1]]
+    iter <- as.integer(unlist(regmatches(iter, gregexpr("[[:digit:]]+", iter))))
+    tmp <- read_delim(paste(mypath, ff, sep="/"), delim="\t", col_names = F) %>% 
+      dplyr::select(-X1, -X3, -X5, -X9, -X11, -X13) %>% 
+      rename(simulated=X2, trait=X4, envir=X6, stat.type=X7, stat.obs=X8, pvalue=X10, SES=X12, exp.med=X14 ) %>% 
+      mutate(stat.type=gsub(pattern = "^r\\(", replacement="", stat.type)) %>%
+      mutate(stat.type=gsub(pattern = "\\)\\=$", replacement="", stat.type)) %>% 
+      mutate(stat.type=gsub(pattern = "_", replacement="\\.", stat.type)) %>% 
+      mutate(trait=gsub(pattern="[[:space:]]+$", replacement="", trait)) %>% 
+      mutate(envir=gsub(pattern="[[:space:]]+$", replacement="", envir)) %>% 
+      mutate(main=iter[[1]]) %>% 
+      mutate(inter=iter[[2]]) %>% 
+      mutate(corr=iter[[3]]) %>% 
+      dplyr::select(main:corr, everything())
+    output <- bind_rows(output, tmp)
+  }
+  
+  outp.summary <- output %>% 
+    dplyr::filter(!stat.type %in% c("TY", "XY", "XY.T", "XY.TR")) %>% 
+    group_by(main, inter, corr, trait, envir, stat.type) %>% 
+    summarize(stat.obs.med=median(stat.obs),
+              power=mean(pvalue<=0.05),
+              SES.med=median(SES), 
+              exp.med.med=median(exp.med), 
+              nsim=n()) %>% 
+    bind_rows(output %>% 
+                dplyr::filter(stat.type %in% c("TY", "XY", "XY.T", "XY.TR")) %>% 
+                group_by(main, inter, corr, trait, stat.type) %>% 
+                summarize(stat.obs.med=median(stat.obs),
+                          power=mean(pvalue<=0.05),
+                          SES.med=median(SES), 
+                          exp.med.med=median(exp.med), 
+                          nsim=n())) %>% 
+    dplyr::select(main:stat.type, nsim, stat.obs.med:exp.med.med) %>% 
+    arrange(stat.type, main, inter, corr, trait, envir)
+  return(outp.summary)
+}
+
+
+get.ntraits <- function(x){ 
+  tmp <- str_split(x, pattern = " ")[[1]]
+  return(length(tmp))
+}
+
+
+#### FIGURE 2 #####
+mypath <- "_data/Experiment_02Mar2020_FactorInteraction&TraitCorr"
+myfiles <- list.files(path=mypath, pattern = "Summary.txt", recursive = T)
+outp.summary <- FormatData(myfiles) 
+
+
+#e41a1c ##red
+#377eb8 ##blue
+#ffffb2 ## yellow
+#4daf4a ##green
+#984ea3 ##violet
+#ff7f00 ##orange
+
+
+mypalette <- palette(c("#e41a1c",  #1 - red)
+        "#ff7f00", #12 - orange
+        "#984ea3", #13 - violet
+        "##ffed6f", #2 - yellow
+        "#4daf4a", #23 - green
+        "#377eb8")) #3 - blue
+
+ggplot(data=outp.summary %>% 
+         mutate(main=main/100) %>% 
+         mutate(corr=factor(corr/10, levels=c(0, 0.4, 0.8), labels=paste0("Correlation = ", c(0, 0.4, 0.8)))) %>% 
+         mutate(inter=factor(inter/10, levels=c(0, 0.3, 0.5), labels=paste0("Interaction = ", c(0, 0.3, 0.5)))) %>% 
+         ungroup() %>% 
+         dplyr::filter(stat.type=="XY") %>% 
+         dplyr::filter(trait %in% c("1", "2", "1 2", "3", "1 3", "2 3")) %>% 
+         mutate(trait=factor(trait, levels=c("1", "2", "3", "1 2", "1 3", "2 3"), 
+                             labels=c("t1", "t2", "tn", "t1 t2", "t1 tn", "t2 tn")))) + #%>%
+         #dplyr::filter(trait %in% c("1", "2", "1 2", "3")) %>%
+         #mutate(inter=as.factor(inter))) + 
+  geom_line(aes(x=main, y=power, group=trait, col=trait)) + 
+  #scale_colour_brewer(palette = "Dark2") + 
+  scale_color_manual("Trait\ncomb.",
+  values=c("#e41a1c",  #1 - red)
+          "#e6ab02", #2 - yellow
+          "#377eb8", #3 - blue
+          "#d95f02", #12 - orange
+          "#984ea3", #13 - violet
+          "#4daf4a" #23 - green
+          ) 
+                      ) +
+  facet_grid(corr~inter) + 
+  theme_bw() + 
+  scale_x_continuous(name="Effect of factor e1 -> trait t1") + 
+  scale_y_continuous(name="Power") + 
+  theme(panel.grid = element_blank())
+
+ggsave(filename="_pics/Fig2_CorrInte_02March.png", width=6, height=5, device="png", dpi = 300, last_plot())
+
+
+
+
+#### FIGURE 3 #####
+mypath <- "_data/Experiment_04Mar2020_TraitNumber"
+myfiles <- list.files(path=mypath, pattern = "Summary.txt", recursive = T)
+outp.summary <- FormatData(myfiles) %>% 
+  rename(ntraits=inter)
+
+## plotting power for XY with corr
+add.t.label <- function(x) {
+  x <- gsub(pattern=" ", replacement=" t", x=x, perl=T)
+  x <- paste0("t", x)
+  return(x)
+}
+
+outp.summary2 <- outp.summary %>% 
+  ungroup() %>% 
+  mutate(main=main/10) %>% 
+  #mutate(corr=factor(corr/10, levels=c(0, 0.4, 0.8), labels=paste0("Correlation = ", c(0, 0.4, 0.8)))) %>% 
+  mutate(ntraits.lab=factor(ntraits, levels=1:5, labels=paste0("n. traits = ", 1:5))) %>%
+  rowwise() %>% 
+  mutate(sel.ntraits=factor(get.ntraits(trait))) %>%  
+  mutate(sel.ntraits.lab=factor(sel.ntraits, levels=levels(sel.ntraits), labels=paste0("Comb. tier = ", levels(sel.ntraits)))) %>%
+  ungroup() %>% 
+  dplyr::filter(stat.type=="XY") %>% 
+  mutate(trait=factor(trait)) %>% 
+  mutate(trait=factor(trait, levels=levels(trait), 
+                      labels=add.t.label(levels(trait))))
+
+### rename null trait to "tn"
+outp.summary2 <- outp.summary2 %>% 
+  mutate(tn.name=paste0("t", ntraits)) %>% 
+  mutate(trait=str_replace(trait, pattern=tn.name, replacement="tn"))
+
+#reorder factors
+outp.summary2$trait <- factor(outp.summary2$trait, levels= c('t1', 't2','t3', 't4', 'tn',
+                                                             't1 t2', 't1 t3','t1 t4','t1 tn','t2 t3','t2 t4','t2 tn','t3 t4','t3 tn', 't4 tn',
+                                                             't1 t2 t3', 't1 t2 t4','t1 t2 tn','t1 t3 t4', 't1 t3 tn','t1 t4 tn','t2 t3 t4','t2 t3 tn', 't2 t4 tn','t3 t4 tn',
+                                                             't1 t2 t3 t4', 't1 t2 t3 tn','t1 t2 t4 tn','t1 t3 t4 tn','t2 t3 t4 tn',
+                                                             't1 t2 t3 t4 tn'))
+
+
+
+hugepalette0 <- c(RColorBrewer::brewer.pal(4, "Dark2"),
+                 gray(0.2),
+                 RColorBrewer::brewer.pal(10, "Paired"),
+                 RColorBrewer::brewer.pal(10, "Set3"),
+                 RColorBrewer::brewer.pal(5, "Pastel1"), "brown")
+#change tone of yellow of t1-t2-t4
+hugepalette0[17] <- "#ccebc5"
+
+hugepalette <- data.frame(trait=c('t1', 't2','t3', 't4', 'tn',
+                                  't1 t2', 't1 t3','t1 t4','t1 tn','t2 t3','t2 t4','t2 tn','t3 t4','t3 tn', 't4 tn',
+                                  't1 t2 t3', 't1 t2 t4','t1 t2 tn','t1 t3 t4', 't1 t3 tn','t1 t4 tn','t2 t3 t4','t2 t3 tn', 't2 t4 tn','t3 t4 tn',
+                                  't1 t2 t3 t4', 't1 t2 t3 tn','t1 t2 t4 tn','t1 t3 t4 tn','t2 t3 t4 tn',
+                                  't1 t2 t3 t4 tn'), 
+                          trait.col=factor(hugepalette0, levels=hugepalette0)) 
+
+outp.summary2 <- outp.summary2 %>% 
+  left_join(hugepalette, by="trait")
+
+fig3 <- ggplot(data=outp.summary2) + 
+  geom_line(aes(x=main, y=power, group=trait, col=trait.col)) + 
+  scale_x_continuous(name="Effect of factor e1 -> trait t1", n.breaks = 4) + 
+  scale_y_continuous(name="Power") + 
+  facet_grid(sel.ntraits.lab~ntraits.lab) +
+  scale_color_identity(guide = "legend", 
+                       labels= hugepalette$trait) +
+  theme_bw() + 
+  theme(panel.grid = element_blank())
+
+
+
+ncols <- c(2,4,4,3,1)
+    
+
+leg.list <- list()
+col.used <- 1
+for(tier in 1:5){
+  outp.summary.tier <- outp.summary2 %>% filter(sel.ntraits==tier)
+  ncombinations <- length(levels(factor(outp.summary.tier$trait)))
+  leg.list[[tier+1]] <- cowplot::get_legend(fig3 %+% outp.summary.tier + 
+                                            guides(col=guide_legend(ncol=ncols[tier], byrow=TRUE))+
+                                            scale_color_identity(name=ifelse(tier==1, "Trait combination - tier 1",paste0("tier - ", tier)), 
+                                                                 guide = "legend", 
+                                                                 labels= hugepalette$trait[col.used:(col.used+ncombinations-1)])
+                                          )
+col.used <- col.used + ncombinations 
+}
+leg.list[[tier+2]] <- NULL
+
+fig3.panel <- cowplot::plot_grid(fig3 + theme(legend.position = "none"), 
+                   cowplot::plot_grid(plotlist = leg.list, nrow = 7, 
+                                      rel_heights = c(0.05, .2,.2,.2,.2,.2, 0.15), align="hv"), 
+                   nrow=1, rel_widths = c(0.6,.4))
+
+
+ggsave(filename="_pics/Fig3_TraitNumber.png", 
+       width=9, height=7, device="png", dpi = 300, fig3.panel)
+
+
+
+#### FIGURE S4 #####
+mypath <- "_data/Experiment_08Jul2020_TraitSuppression&TraitCorr"
+myfiles <- list.files(path=mypath, pattern = "Summary.txt", recursive = T)
+outp.summary <- FormatData(myfiles) 
+
+
+ggplot(data=outp.summary %>% 
+         mutate(main=main/100) %>% 
+         mutate(corr=factor(corr/10, levels=c(0, 0.4, 0.8), labels=paste0("Correlation = ", c(0, 0.4, 0.8)))) %>% 
+         mutate(inter=factor(inter/10, levels=c(0, 0.3, 0.5), labels=paste0("Interaction = ", c(0, 0.3, 0.5)))) %>% 
+         ungroup() %>% 
+         dplyr::filter(stat.type=="XY") %>% 
+         dplyr::filter(trait %in% c("1", "2", "3", "4")) %>% 
+         mutate(trait=factor(trait, levels=c("1", "2", "3", "4"), 
+                             labels=c("t1", "t2", "t3","tn")))) + #%>%
+  #dplyr::filter(trait %in% c("1", "2", "1 2", "3")) %>%
+  #mutate(inter=as.factor(inter))) + 
+  geom_line(aes(x=main, y=power, group=trait, col=trait)) + 
+  #scale_colour_brewer(palette = "Dark2") + 
+  #scale_color_manual("Trait\ncomb.",
+  scale_color_manual("Trait",
+                     values=c("#e41a1c",  #1 - red)
+                              "#e6ab02", #2 - yellow
+                              "#4daf4a", #23 - green
+                              "#377eb8", #3 - blue
+                              "#d95f02", #12 - orange
+                              "#984ea3" #13 - violet
+                     ) 
+  ) +
+  facet_grid(corr~inter) + 
+  theme_bw() + 
+  scale_x_continuous(name="Effect of factor e1 -> trait t1") + 
+  scale_y_continuous(name="Power") + 
+  theme(panel.grid = element_blank())
+
+ggsave(filename="_pics/FigS4_Extra_CorrInte_08Jul20.png", width=6, height=5, device="png", dpi = 300, last_plot())
+
+
+
+
+
+
+
+
+
+
+
+
+
-- 
GitLab