03_Extra_Figures_Simulations.R

## 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 S4 #####
mypath <- "_data/Experiment_08Jul2020_TraitSuppression&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", "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/Fig2_Extra_CorrInte_08Jul20.png", width=6, height=5, device="png", dpi = 300, last_plot())