# =============================== Load Library ===============================
library('R.matlab')
library('reshape2')
library('doBy')
library('grid')
library('ggplot2')
library('dplyr')
library('tidyr') 
library('caret')
library('e1071')

# ------------------------------ Read real data ------------------------------
timeGen = readRDS('timeGen_s4.RDS')

# ----------------------------- Process real data ----------------------------
timeGen_mean = timeGen %>% 
  group_by(timeGen, time, TestTime) %>% 
  summarise(acc = mean(acc)) %>% 
  ungroup()

# --------------------- Make Grid reference to read data ---------------------
iGrid = expand.grid(ic = 1:36, ip = 1:100) %>% as.tbl()

condiPair = data.frame(c1 = c(111, 115, 112, 116, 114, 118),
                       c2 = c(125, 121, 126, 122, 128, 124))
condiPair_h = condiPair

geneTable = condiPair_h %>% 
  group_by(c1, c2) %>% 
  do(condiPair_h %>% rename(g1 = c1, g2= c2)) %>%
  ungroup()

# ----------------------------------------------------------------------------
geneTable$ind = 1:nrow(geneTable)
iPlotCondi = geneTable %>% filter(c1 == condiPair$c1[6], c2 == condiPair$c2[6])
iPlotCondi = iPlotCondi$ind



# ------------------------------ Read perm data ------------------------------
timeGen.perm = iGrid %>% 
  filter(ic %in% iPlotCondi) %>%
  rowwise() %>% 
  do( readRDS(sprintf('perm\\timeGen_perm_c%d_r%d.RDS', .$ic, .$ip)) )


# ---------------------------- calculate perm data ---------------------------
timeGen.perm.summary = timeGen.perm %>%
  ungroup %>%
  group_by(time, TestTime, timeGen) %>%
  summarise(CI = quantile(acc, .95), sd = sd(acc), acc_perm = mean(acc))


# ----------------------------------------------------------------------------
timeGen.plot = timeGen.perm.summary %>% 
  left_join(timeGen_mean) %>% 
  ungroup() %>%
  mutate(isSig = as.integer(acc > CI))

# ----------------------------------------------------------------------------
ggplot(timeGen.plot, aes(TestTime, time, fill = acc)) + 
  geom_tile( colour = "black" ) + 
  # scale_fill_gradientn(colours = colorRamps::matlab.like2(100)) + 
  geom_tile( aes(fill=isSig)) +
  facet_wrap(~timeGen, ncol = 3) + 
  coord_fixed()