### R code from vignette source 'vignettes/HiTC/inst/doc/HiTC.Rnw'

###################################################
### code chunk number 1: head
###################################################
library(HiTC)
showClass("HTCexp")


###################################################
### code chunk number 2: import5c
###################################################
exDir <- system.file("extdata", package="HiTC")
GM12878<-import.my5C(file.path(exDir,"nsmb.1936-S5.txt"), 
                 xgi.bed=file.path(exDir,"Bau_GM12878_REV.bed"), ygi.bed=file.path(exDir,"Bau_GM12878_FOR.bed")) 

## List of HTCexp objects describing the 5C dataset                                                                 
detail(GM12878)


###################################################
### code chunk number 3: qc5c (eval = FALSE)
###################################################
## CQC(GM12878)


###################################################
### code chunk number 4: QCeval
###################################################
png(file="HiTC-qc5c.png", res=300, units="in", width=5, height=5)
CQC(GM12878)
graphics.off()


###################################################
### code chunk number 5: mapC5c (eval = FALSE)
###################################################
## mapC(GM12878$chr16chr16)


###################################################
### code chunk number 6: norm5CEval
###################################################
png(file="HiTC-mapC.png", res=300, units="in", width=5, height=5)
mapC(GM12878$chr16chr16)
graphics.off()


###################################################
### code chunk number 7: bin5C
###################################################

## Binning of 5C interaction map
GM12878.binned <- binningC(GM12878$chr16chr16, binsize=50000, step=5)
mapC(GM12878.binned, show.na=TRUE)


###################################################
### code chunk number 8: norm5Cexp (eval = FALSE)
###################################################
## ## Look at exptected counts
## GM12878exp <- getExpectedCounts(GM12878$chr16chr16, stdev=TRUE, plot=TRUE, bin=0.001, span=0.3)


###################################################
### code chunk number 9: norm5CexpEval
###################################################
png(file="HiTC-expint.png", res=300, units="in", width=6, height=6)
GM12878exp <- getExpectedCounts(GM12878$chr16chr16, bin=0.001, span=0.3, stdev=TRUE, plot=TRUE)
graphics.off()


###################################################
### code chunk number 10: norm5Cznorm
###################################################
GM12878z <- normPerZscore(GM12878$chr16chr16, span=0.3)
mapC(GM12878z)


###################################################
### code chunk number 11: annot5C
###################################################
Refgene <- readBED(file.path(exDir,"refseq_hg19_chr16_1_500000.bed"))
mapC(GM12878.binned, 
     giblocs=list(RefSeqGene=Refgene$Refseq_Gene), view=2)


###################################################
### code chunk number 12: comp5C
###################################################
K562<-import.my5C(file.path(exDir,"nsmb.1936-S6.txt"), 
                 xgi.bed=file.path(exDir,"Bau_K562_REV.bed"), ygi.bed=file.path(exDir,"Bau_K562_FOR.bed")) 
K562.binned <- binningC(K562$chr16chr16, binsize=50000, step=5)

mapC(GM12878.binned, K562.binned,
     giblocs=list(RefSeqGene=Refgene$Refseq_Gene))


###################################################
### code chunk number 13: importhic
###################################################
## Load Dekker et al. Chromosome 14 data (from GEO GSE18199)
hiC<-import.my5C(file.path(exDir,"HIC_gm06690_chr14_chr14_1000000_obs.txt"), 
                 xgi.bed=file.path(exDir,"GSE18199_gm06690_chr14_chr14_1Mb.bed")) 
hiC <- hiC$chr14chr14
detail(hiC)


###################################################
### code chunk number 14: mapChic
###################################################
## Extract region of interest and plot the interaction map
hiC <- extractRegion(hiC,from=1.8e+07, to=106368584)
mapC(hiC, maxrange=100)


###################################################
### code chunk number 15: mapNormhic
###################################################
## Data Normalisation by Expected number of Counts
hiCnorm <- normPerExpected(hiC)
mapC(hiCnorm, log.data=TRUE)


###################################################
### code chunk number 16: mapCorhic
###################################################
## Correlation Map of Chromosome 14
mapC(cor(intdata(hiCnorm)), maxrange=1, minrange=-1, 
     col.pos=c("black", NA, "red"), col.neg=c("black",NA, "blue"))


###################################################
### code chunk number 17: sessionInfo
###################################################
toLatex(sessionInfo(), locale=FALSE)