library(ggcyto)
## Loading required package: grid
## Loading required package: ggplot2
## Loading required package: flowCore
## Loading required package: ncdfFlow
## Loading required package: flowViz
## Loading required package: lattice
## Loading required package: RcppArmadillo
## Loading required package: BH
## Loading required package: flowWorkspace
## Loading required package: gridExtra
## Warning: replacing previous import by 'ncdfFlow::rbind2' when loading
## 'ggcyto'
library(flowIncubator)
## Loading required package: QUALIFIER
## Loading required package: data.table
## Loading required package: reshape
path <- "~/rglab/workspace/flowWorkspace/wsTestSuite/gatingML/"
# xmlfile <- file.path(path, "Merck/CytExp_10623_Gates_v5.xml")
xmlfile <- system.file("extdata/cytotrol_tcell_cytobank.xml", package = "flowIncubator")

g <- read.gatingML.cytobank(xmlfile)
g
## --- Gating hieararchy parsed from GatingML: 
##  with  5  populations defined
#tree accessors
getNodes(g)
## GateSet_686424 GateSet_686425 GateSet_686426 GateSet_686427 GateSet_686428 
##   "not debris"     "singlets"          "CD3"          "CD8"          "CD4"
getParent(g, "GateSet_686425")
## [1] "GateSet_686424"
getChildren(g, "GateSet_686425")
## [1] "GateSet_686426"
plot(g)

#read FCS
# fcsFiles <- list.files(pattern = "\\.fcs", file.path(path, "Merck"), full = T)
fcsFiles <- list.files(pattern = "CytoTrol", system.file("extdata", package = "flowWorkspaceData"), full = T)
fs_raw <- read.ncdfFlowSet(fcsFiles)
## All FCS files have the same following channels:
## FSC-A
## FSC-H
## FSC-W
## SSC-A
## B710-A
## R660-A
## R780-A
## V450-A
## V545-A
## G560-A
## G780-A
## Time
## write CytoTrol_CytoTrol_1.fcs to empty cdf slot...
## write CytoTrol_CytoTrol_2.fcs to empty cdf slot...
## done!
# save_ncfs(fs_raw, path = "/loc/no-backup/mike/fs_raw")

# fs_raw <- load_ncfs("/loc/no-backup/mike/fs_raw")
# fs <- clone.ncdfFlowSet(fs_raw,isEmpty = F)

#compensation
comp <- getCompensationMatrices(g)
fs <- compensate(fs_raw, comp)
## copying data slice:CytoTrol_CytoTrol_1.fcs
## write CytoTrol_CytoTrol_1.fcs to empty cdf slot...
## copying data slice:CytoTrol_CytoTrol_2.fcs
## write CytoTrol_CytoTrol_2.fcs to empty cdf slot...
## updating CytoTrol_CytoTrol_1.fcs...
## updating CytoTrol_CytoTrol_2.fcs...
#transforamation
trans <- getTransformations(g)
fs <- transform(fs, trans)
## copying data slice:CytoTrol_CytoTrol_1.fcs
## write CytoTrol_CytoTrol_1.fcs to empty cdf slot...
## copying data slice:CytoTrol_CytoTrol_2.fcs
## write CytoTrol_CytoTrol_2.fcs to empty cdf slot...
## updating CytoTrol_CytoTrol_1.fcs...
## updating CytoTrol_CytoTrol_2.fcs...
#check the trans and comps outcome
ggcyto(fs, aes(x = "CD4"), subset = "root") + geom_density()

ggcyto(fs, aes(x = "CD4", y = "CD8"), subset = "root") + geom_hex()

#apply gates to gatingset
gs <- GatingSet(fs)
## ..done!
gating(g, gs)
## not debris
## singlets
## CD3
## CD4
## CD8
## ..done!
#check the gating results
plotGate(gs[[1]])

#load stats from cytobank
cytobank_counts <- read.csv(system.file("extdata/cytotrol_tcell_cytobank_counts.csv", package = "flowIncubator"), skip = 7, stringsAsFactors = F)
#convert it to the same format
colnames(cytobank_counts)[1] <- "Population"
cytobank_counts <- plyr::arrange(cytobank_counts, Population)

#load openCyto stats
openCyto_counts <- getPopStats(gs, statType = "count")
openCyto_counts <- reshape2::dcast(openCyto_counts, Population ~ name, value.var = "Count")

#compare two
all.equal(cytobank_counts, openCyto_counts)
## [1] TRUE