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

###################################################
### code chunk number 1: load
###################################################
library(ggbio)
data(hg19IdeogramCyto, package = "biovizBase")
data(hg19Ideogram, package = "biovizBase")
library(GenomicRanges)


###################################################
### code chunk number 2: simul_gr
###################################################
library(biovizBase)
gr <- GRanges(rep(c("chr1", "chr2"), each = 5),
              IRanges(start = rep(seq(1, 100, length = 5), times = 2),
                      width = 50))

autoplot(gr, aes(fill = seqnames))


###################################################
### code chunk number 3: coord-genome
###################################################
autoplot(gr, coord = "genome", aes(fill = seqnames))


###################################################
### code chunk number 4: is
###################################################
gr.t <- transformToGenome(gr)
head(gr.t)
is_coord_genome(gr.t)
metadata(gr.t)$coord


###################################################
### code chunk number 5: simul_snp
###################################################
chrs <- as.character(levels(seqnames(hg19IdeogramCyto)))
seqlths <- seqlengths(hg19Ideogram)[chrs]
set.seed(1)
nchr <- length(chrs)
nsnps <- 100
gr.snp <- GRanges(rep(chrs,each=nsnps),
                  IRanges(start =
                          do.call(c, lapply(chrs, function(chr){
                            N <- seqlths[chr]
                            runif(nsnps,1,N)
                          })), width = 1),
                  SNP=sapply(1:(nchr*nsnps), function(x) paste("rs",x,sep='')),
                  pvalue =  -log10(runif(nchr*nsnps)),
                  group = sample(c("Normal", "Tumor"), size = nchr*nsnps,
                    replace = TRUE)
                  )
genome(gr.snp) <- "hg19"
gr.snp


###################################################
### code chunk number 6: shorter
###################################################
seqlengths(gr.snp)
nms <- seqnames(seqinfo(gr.snp))
nms.new <- gsub("chr", "", nms)
names(nms.new) <- nms
gr.snp <- renameSeqlevels(gr.snp, nms.new)
seqlengths(gr.snp)


###################################################
### code chunk number 7: unorder
###################################################
autoplot(gr.snp, coord = "genome", geom = "point", aes(y = pvalue), space.skip = 0.01)


###################################################
### code chunk number 8: sort
###################################################
gr.snp <- keepSeqlevels(gr.snp, c(1:22, "X", "Y"))
autoplot(gr.snp, coord = "genome", geom = "point", aes(y = pvalue), space.skip = 0.01)


###################################################
### code chunk number 9: with_seql
###################################################
names(seqlths) <- gsub("chr", "", names(seqlths))
seqlengths(gr.snp) <- seqlths[names(seqlengths(gr.snp))]
autoplot(gr.snp, coord = "genome", geom = "point", aes(y = pvalue), space.skip = 0.01)


###################################################
### code chunk number 10: line
###################################################
autoplot(gr.snp, coord = "genome", geom = "line", aes(y = pvalue, group = seqnames,
                                     color = seqnames))


###################################################
### code chunk number 11: plotGrandLinear
###################################################
plotGrandLinear(gr.snp, aes(y = pvalue))


###################################################
### code chunk number 12: morecolor1
###################################################
plotGrandLinear(gr.snp, aes(y = pvalue, color = seqnames))


###################################################
### code chunk number 13: morecolor2
###################################################
plotGrandLinear(gr.snp, aes(y = pvalue), color = c("green", "deepskyblue"))


###################################################
### code chunk number 14: morecolor3
###################################################
plotGrandLinear(gr.snp, aes(y = pvalue), color = c("green", "deepskyblue", "red"))


###################################################
### code chunk number 15: morecolor4
###################################################
plotGrandLinear(gr.snp, aes(y = pvalue), color = "red")


###################################################
### code chunk number 16: cutoff
###################################################
plotGrandLinear(gr.snp, aes(y = pvalue), cutoff = 3, cutoff.color = "blue", cutoff.size = 4)


###################################################
### code chunk number 17: cutoff-low (eval = FALSE)
###################################################
## plotGrandLinear(gr.snp, aes(y = pvalue)) + geom_hline(yintercept = 3, color = "blue", size = 4)


###################################################
### code chunk number 18: longer
###################################################
## let's make a long name
nms <- seqnames(seqinfo(gr.snp))
nms.new <- paste("chr00000", nms, sep = "")
names(nms.new) <- nms
gr.snp <- renameSeqlevels(gr.snp, nms.new)
seqlengths(gr.snp)


###################################################
### code chunk number 19: rotate
###################################################
plotGrandLinear(gr.snp, aes(y = pvalue)) + opts(axis.text.x=theme_text(angle=-90, hjust=0))


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