### R code from vignette source 'CNEr.Rnw'

###################################################
### code chunk number 1: initialize
###################################################
options(width=70)


###################################################
### code chunk number 2: rtl-Axt
###################################################
library(CNEr)
axtFilesHg19DanRer7 <- file.path(system.file("extdata", package="CNEr"), 
                                "hg19.danRer7.net.axt")
axtFilesDanRer7Hg19 <- file.path(system.file("extdata", package="CNEr"), 
                                "danRer7.hg19.net.axt")


###################################################
### code chunk number 3: rtl-Axt (eval = FALSE)
###################################################
## axtHg19DanRer7 <- readAxt(axtFilesHg19DanRer7)
## axtDanRer7Hg19 <- readAxt(axtFilesDanRer7Hg19)


###################################################
### code chunk number 4: <loadAxt
###################################################
data(axtHg19DanRer7)
data(axtDanRer7Hg19)


###################################################
### code chunk number 5: <<loadAxt
###################################################
axtHg19DanRer7
axtDanRer7Hg19


###################################################
### code chunk number 6: rtl-Bed
###################################################
bedHg19Fn <- file.path(system.file("extdata", package="CNEr"), 
                      "filter_regions.hg19.bed")
bedHg19 <- readBed(bedHg19Fn)
bedHg19
bedDanRer7Fn <- file.path(system.file("extdata", package="CNEr"), 
                         "filter_regions.danRer7.bed")
bedDanRer7 <- readBed(bedDanRer7Fn)
bedDanRer7


###################################################
### code chunk number 7: rtl-chromSizes (eval = FALSE)
###################################################
## qSizesHg19 <- fetchChromSizes("hg19")
## qSizesDanRer7 <- fetchChromSizes("danRer7")


###################################################
### code chunk number 8: chromSizesData
###################################################
data(qSizesHg19)
data(qSizesDanRer7)


###################################################
### code chunk number 9: showchromSizesData
###################################################
qSizesHg19
qSizesDanRer7


###################################################
### code chunk number 10: rtl-CNEScan
###################################################
## axt, GRanges objects as input
CNEHg19DanRer7 <- ceScan(axts=axtHg19DanRer7, tFilter=bedHg19,
                         qFilter=bedDanRer7, qSizes=qSizesDanRer7,
                         thresholds=c("45_50", "48_50", "49_50"))
CNEDanRer7Hg19 <- ceScan(axts=axtDanRer7Hg19, tFilter=bedDanRer7,
                         qFilter=bedHg19, qSizes=qSizesHg19,
                         thresholds=c("45_50", "48_50", "49_50"))

## axt and bed files as input
CNEHg19DanRer7 <- ceScan(axts=axtFilesHg19DanRer7, tFilter=bedHg19Fn,
                         qFilter=bedDanRer7Fn, qSizes=qSizesDanRer7,
                         thresholds=c("45_50", "48_50", "49_50"))
CNEDanRer7Hg19 <- ceScan(axts=axtFilesDanRer7Hg19, tFilter=bedDanRer7Fn,
                         qFilter=bedHg19Fn, qSizes=qSizesHg19,
                         thresholds=c("45_50", "48_50", "49_50"))


###################################################
### code chunk number 11: rtl-CNEScan
###################################################
#data(CNEHg19DanRer7)
lapply(CNEHg19DanRer7, head)


###################################################
### code chunk number 12: rtl-CNEMerge
###################################################
#data(CNEDanRer7Hg19)
cneMergedDanRer7Hg19 <- mapply(cneMerge, CNEDanRer7Hg19, CNEHg19DanRer7, 
                               SIMPLIFY=FALSE)
lapply(cneMergedDanRer7Hg19, head)


###################################################
### code chunk number 13: rtl-CNERealignment (eval = FALSE)
###################################################
## assemblyHg19Twobit <- "/Users/gtan/CSC/CNEr/2bit/hg19.2bit"
## assemblyDanRer7Twobit <- "/Users/gtan/CSC/CNEr/2bit/danRer7.2bit"
## cneBlatedDanRer7Hg19 <- list()
## for(i in 1:length(cneMergedDanRer7Hg19)){
##   cneBlatedDanRer7Hg19[[names(cneMergedDanRer7Hg19)[i]]] <-
##     blatCNE(cneMergedDanRer7Hg19[[i]], 
##             as.integer(sub(".+_.+_\\d+_", "", names(cneMergedDanRer7Hg19)[i])), 
##             cutoffs1=8L, cutoffs2=8L, 
##             assembly1Twobit=assemblyDanRer7Twobit, 
##             assembly2Twobit=assemblyHg19Twobit,
##             blatBinary="blat")
## }


###################################################
### code chunk number 14: ceScanOneStep (eval = FALSE)
###################################################
## assemblyHg19Twobit <- "/Users/gtan/CSC/CNEr/2bit/hg19.2bit"
## assemblyDanRer7Twobit <- "/Users/gtan/CSC/CNEr/2bit/danRer7.2bit"
## finalCNE <- ceScanOneStep(axt1=axtHg19DanRer7, filter1=bedHg19, 
##                           sizes1=qSizesHg19, assembly1="hg19", 
##                           twoBit1=assemblyHg19Twobit,
##                           axt2=axtDanRer7Hg19, filter2=bedDanRer7,
##                           sizes2=qSizesDanRer7, assembly2="danRer7",
##                           twoBit2=assemblyDanRer7Twobit,
##                           thresholds=c("45_50", "48_50", "49_50"),
##                           blatBinary="blat", blatCutoff1=8L, blatCutoff2=8L)


###################################################
### code chunk number 15: saveCNE
###################################################
## on individual tables
dbName <- tempfile()
data(cneBlatedDanRer7Hg19)
for(i in 1:length(cneBlatedDanRer7Hg19)){
  tableName <- paste("danRer7_hg19", names(cneBlatedDanRer7Hg19)[i],
                      sep="_")
	saveCNEToSQLite(cneBlatedDanRer7Hg19[[i]], dbName, tableName, overwrite=TRUE)
}

## on CNE class
data(finalCNE)
saveCNEToSQLite(finalCNE, dbName=dbName, overwrite=TRUE)


###################################################
### code chunk number 16: queryCNE
###################################################
chr <- "chr11"
start <- 31000000L
end <-  33000000L
minLength <- 50L
tableName <- "danRer7_hg19_45_50"
fetchedCNERanges <- readCNERangesFromSQLite(dbName, tableName, chr, 
                                           start, end, whichAssembly="L",
                                           minLength=minLength)


###################################################
### code chunk number 17: queryUCSC (eval = FALSE)
###################################################
## library(Gviz)
## genome <- "hg19"
## chr <- "chr11"
## start <- 31000000L
## end <-  33000000L
## axisTrack <- GenomeAxisTrack()
## ideoTrack <- IdeogramTrack(genome=genome, chromosome=chr)
## cpgIslands <- UcscTrack(genome=genome, chromosome=chr,
##                        track="cpgIslandExt", from=start, to=end,
##                        trackType="AnnotationTrack", start="chromStart",
##                        end="chromEnd", id="name", shape="box",
##                        fill="#006400", name="CpG Islands")
## refGenes <- UcscTrack(genome="hg19", chromosome=chr,
##                      track="refGene", from=start, to=end,
##                      trackType="GeneRegionTrack", rstarts="exonStarts",
##                      rends="exonEnds", gene="name2", symbol="name2",
##                      transcript="name", strand="strand", fill="#8282d2",
##                      name="refSeq Genes", collapseTranscripts=TRUE)
## biomTrack <- BiomartGeneRegionTrack(genome="hg19", chromosome=chr, 
##                                    start=start , end=end, name="Ensembl")


###################################################
### code chunk number 18: loadAnnotation
###################################################
data(axisTrack)
data(ideoTrack)
data(cpgIslands)
data(refGenes)


###################################################
### code chunk number 19: plotAnnotation
###################################################
library(Gviz)
plotTracks(list(axisTrack, ideoTrack, cpgIslands, refGenes), 
           collapseTranscripts=TRUE, shape="arrow",
           showId=TRUE, transcriptAnnotation="symbol")


###################################################
### code chunk number 20: CNEDensity
###################################################
dbName <- file.path(system.file("extdata", package="CNEr"),
                    "cne.sqlite")
windowSize <- 300L
minLength <- 50L
cneHg19DanRer7_45_50 <- 
  CNEDensity(dbName=dbName, 
             tableName="danRer7_hg19_45_50", 
             assembly1="hg19", chr=chr, start=start,
             end=end, windowSize=windowSize, 
             minLength=minLength)
cneHg19DanRer7_48_50 <-
  CNEDensity(dbName=dbName, 
             tableName="danRer7_hg19_48_50", 
             assembly1="hg19", chr=chr, start=start,
             end=end, windowSize=windowSize, 
             minLength=minLength)
cneHg19DanRer7_49_50 <-
  CNEDensity(dbName=dbName, 
             tableName="danRer7_hg19_49_50", 
             assembly1="hg19", chr=chr, start=start,
             end=end, windowSize=windowSize, 
             minLength=minLength)


###################################################
### code chunk number 21: GvizDataTrack
###################################################
data(cneHg19DanRer7_45_50)
data(cneHg19DanRer7_48_50)
data(cneHg19DanRer7_49_50)
genome <- "hg19"
chr <- "chr11"
start <- 31000000L
end <-  33000000L
#axisTrack = GenomeAxisTrack()
#ideoTrack = IdeogramTrack(genome=genome, chromosome=chr)
strand <- "+"
dataMatrix <- cneHg19DanRer7_45_50
dTrack1 <- DataTrack(start=dataMatrix[ ,1], end=dataMatrix[ ,1], 
                     data=dataMatrix[ ,2], chromosome=chr, strand=strand, 
                     genome=genome, type="horiz", horizon.scale=0.1, 
                     fill.horizon=c("#B41414", "#E03231", "#F7A99C", 
                                    "yellow", "orange", "red"), 
                     name="danRer7 45/50")
dataMatrix <- cneHg19DanRer7_48_50
dTrack2 <- DataTrack(start=dataMatrix[ ,1], end=dataMatrix[ ,1], 
                     data=dataMatrix[ ,2], chromosome=chr, strand=strand, 
                     genome=genome, type="horiz", horizon.scale=0.1, 
                     fill.horizon=c("#B41414", "#E03231", "#F7A99C", 
                                    "yellow", "orange", "red"), 
                     name="danRer7 48/50")
dataMatrix <- cneHg19DanRer7_49_50
dTrack3 <- DataTrack(start=dataMatrix[ ,1], end=dataMatrix[ ,1], 
                     data=dataMatrix[ ,2], chromosome=chr, strand=strand, 
                     genome=genome, type="horiz", horizon.scale=0.1, 
                     fill.horizon=c("#B41414", "#E03231", "#F7A99C", 
                                    "yellow", "orange", "red"), 
                     name="danRer7 49/50")


###################################################
### code chunk number 22: plotCNE
###################################################
plotTracks(list(axisTrack, ideoTrack, cpgIslands, 
                refGenes, 
                dTrack1, dTrack2, dTrack3),
           collapseTranscripts=TRUE, shape="arrow",
           showId=TRUE, transcriptAnnotation="symbol")


###################################################
### code chunk number 23: session-info
###################################################
  sessionInfo()