# Real data study - Part 1


# Set working directory and path to store results to:
setwd("~/")
savepath <- "~/"


# Accession numbers of datasets to use:

accnumber <- c("E-GEOD-60542",  
 "E-GEOD-71571",  
 "E-GEOD-64415",  
 "E-GEOD-48278",  
 "E-GEOD-58697",  
 "E-GEOD-68720",  
 "E-GEOD-67851",  
 "E-GEOD-64857",  
 "E-GEOD-53224",  
 "E-GEOD-49697",  
 "E-GEOD-28654",  
 "E-GEOD-54543",  
 "E-GEOD-63626",  
 "E-GEOD-57065",  
 "E-GEOD-61804",  
 "E-GEOD-33532",  
 "E-GEOD-19722",  
 "E-GEOD-59312",  
 "E-GEOD-39084",  
 "E-GEOD-60028",  
 "E-GEOD-45670",  
 "E-GEOD-38662",  
 "E-GEOD-58558",  
 "E-GEOD-49243",  
 "E-GEOD-54837",  
 "E-GEOD-50774",  
 "E-GEOD-46474",  
 "E-GEOD-48780",  
 "E-GEOD-48350",  
 "E-GEOD-53890",  
 "E-GEOD-45216",  
 "E-GEOD-53987",  
 "E-GEOD-29623",  
 "E-GEOD-45867")

 
# Make table of settings: 
 
setting <- c("sepnorm", "addonnorm")
batchremmethod <- c("none", "fabatch", "combat", "fsvafast", "fsvaexact", "meancenter", 
  "standardize", "ratioa", "ratiog")
traintestsetsize <- c("trainsmalltestverysmall", "trainsmalltestsmall", 
 "trainsmalltestlarge", "trainlargetestverysmall", "trainlargetestsmall", 
 "trainlargetestlarge")
  
scenariogrid <- expand.grid(accnumbertest=accnumber, accnumbertrain=accnumber, traintestsetsize=traintestsetsize,
  batchremmethod=batchremmethod, setting=setting, stringsAsFactors = FALSE)[,5:1]
  
set.seed(1234)
seeds <- sample(1000:10000000, size=length(accnumber)^2)

scenariogrid$seed <- rep(seeds, times=length(batchremmethod)*length(setting)*length(traintestsetsize))

scenariogrid <- scenariogrid[apply(scenariogrid, 1, function(x) x[4]!=x[5]),]



# Randomly permute rows of the table containing the settings.
# This is performed to ensure a comparable computational burdens of the
# four jobs corresponding to the four sub-scripts which perform
# the study.

set.seed(1234)
reorderind <- sample(1:nrow(scenariogrid))
scenariogrid <- scenariogrid[reorderind,]



# Subset the table containing the settings to contain only
# those to be performed by this script:
 
njobs <- 4

numberits <- floor(nrow(scenariogrid)/njobs)
nrest <- nrow(scenariogrid) - njobs*numberits

nitervec <- rep(numberits, njobs)

if(nrest > 0)
  nitervec[1:nrest] <- nitervec[1:nrest] + 1
 
startinds <- c(1, cumsum(nitervec)+1)
startinds <- startinds[-length(startinds)]
stopinds <- cumsum(nitervec)

indjob <- 1
scenariogrid <- scenariogrid[startinds[indjob]:stopinds[indjob],] 
 

# Source the functions to perform the calculations on the server:

source(paste(savepath, "ICSVCodeData/Functions/FunctionsForClusters.R", sep=""))


# Export all objects to the CPUs:

sfExportAll()


# Evaluate the settings:

Results <- sfLapply(1:nrow(scenariogrid), function(z) try({evaluatesetting(z)}))

# Stop the cluster:

sfStop()


# Save the results:

save(Results, file=paste(savepath, "ICSVCodeData/Results/ResultsLargeDataStudy_1.Rda", sep=""))


# Clear the workspace:

rm(list=ls());gc()
